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I.  INTRODUCTION 


Nonorthogonal  wavelet  transforms  play  an  important  role  in  signal  processing  by  offer¬ 
ing  much  greater  flexibility  in  the  choice  of  waveform  than  their  orthogonal  counterparts. 
In  particular,  they  allow  adjustment  of  the  relative  resolution  in  time  and  scale.  For  ex¬ 
ample,  decreasing  the  parameter  /3  in  the  Morlet  wavelet  i/’(0  =  ^ increases 

frequency  resolution  while  decreasing  time  resolution  (cf.  [l],  [2],  Appendix  A).  We  are 
coi  cerned  with  discrete  transforms,  however,  and  increased  resolution  does  little  good 
unless  it  is  accompanied  by  a  sufficiently  fine  output  grid.  Such  output  requires  an  undec¬ 
imated  transform  and  the  addition  of  voices  (i.e.,  suboctave  sampling),  which  is  generally 
inconsistent  with  orthogonality.^  Moreover,  many  applications  require  this  finer  sampling 
everx  when  the  shape  of  a  relatively  broadband  mother  wavelet  is  adequate. 

The  standard  inversion  procedure  for  discrete  nonorthogonal  transforms  is  a  finite 
expansion  in  terms  of  the  analyzing  wavelet  [1],  [3],  [4].  Formally,  it  is  based  on  the 
theory  of  frames  and  can  be  thought  of  as  a  a  discretization  of  the  corresponding  con¬ 
tinuous  inversion  formula.  From  this  point  of  view,  several  approximations  are  involved: 
the  wavelet  coefficients  themselves,  the  partial  sum,  and  the  use  of  the  analyzing  wavelets 
rather  than  their  duals.  While  a  partial  expansion  works  quite  well  for  many  (sufficiently 
oscillatory)  signals,  it  fa.ils  to  achieve  good  accuracy  or  requires  an  excessive  number  of 
scales  for  others.  Unfortunately,  the  analyses  found  in  the  literature  focus  on  frame  bounds 
rather  than  the  quality  of  finite  discrete  implementations,  and  generally  only  treat  rela¬ 
tively  broadband  wavelets  (i.e.,  /3  about  2/3).  In  short,  there  is  a  clear  need  for  a  more 
careful  examination  of  the  invertibility  of  these  transforms.  This  paper  provides  several 
alternative  algorithms  for  inversion  of  the  discrete  wavelet  transform  and  compares  them 
in  the  crise  of  Morlet  wavelets.  In  the  process,  both  practical  and  theoretical  issues  for  the 
inversion  of  nonorthogonal  discrete  wavelet  transforms  are  discussed. 

Our  ultimate  goal  is  to  formulate  the  problem  in  an  entirely  discrete  context,  unifying 
the  various  inverses  under  a  common  filter  bank  stiucture.  The  continuous  case  shall  serve 
as  guide  and  interpreter.  The  remainder  of  the  introduction  is  devoted  to  elaborating 
this  point  of  view.  Section  II  formally  introduces  our  notation,  and  Section  III  describes 
the  dual  wavelets  in  a  general  setting,  providing  some  results  particular  to  undecimated 
transforms.  Section  IV  discretizes  the  continuous  case.  Section  V  begins  with  a  general 
treatment  of  the  discrete  inverse  wavelet  transform  including  a  discussion  of  appropriate 
filter  constraints.  Then,  three  particularly  relevant  inverses  are  pre.sented  and  analyzed. 
Section  VI  provides  a  numerical  comparison  of.  the  three  inverses  and  the  frame  approxi¬ 
mation. 


^  This  statenfiftit  is  no  longer  true  if  one  admits  bases  generated  by  more  than  one  mother  wavelet. 
Wavelet  packets  also  present  a  means  of  obtaining  higher  resolution.  These  methods  suffer  various  draw¬ 
backs,  however,  which  are  outside  the  scope  of  this  paper. 


Let  \b(t)  be  a  square  integrable  function  satisfying  the  admissibility  condition 


where 


(1.1) 


(1.2) 


is  the  Fourier  transform  of  We  define  the  continuous  wavelet  transform  of  a  .signal 

s(t)  by 

iy(a,  b)  =  ’  (1-3) 

where  the  bar  indicates  complex  conjugate.  For  orthonormal  wavelets,  the  following  ex¬ 
pansion  is  valid: 


sit)  =  ^  Wi2\Tn)^i,Qr  -n^ 

=  V  W..„  V)..„  . 


»,n 


(1.4) 


One  might  hope  for  a  similar  expansion  for  nonorthogonal  wavelets.  Let  us  consider 
discretizing  the  continuous  case.  If  rpit)  is  reed  or  analytic  and  satisfies  (l-l),  then  the 
signal  maj'  be  recovered  by  the  standard  inversion  formula,  [1],  [3]^ 

^ r  ^  ly 

Setting  a  =  2‘  and  b  =  2'n  in  (1.5),  we  obtain  (cf.  Section  IV) 

~  lV(2‘,2*u)V’ ,  (1.6) 

^  i,n 

which  is  essentially  the  same  expansion  as  (1.4).  Unfortunately  (1.6)  is  not,  in  general, 
correct.  The  theory  of  frames  [3] -[4]  tells  us  that 

s{t)  =  Re  Y,  W.,„  i^nit)  ,  (1.7) 

i,n 

■^If  —  0  tor  O'  <  0,  then  the  formula  is  valid  for  complex  ip  provided  one  takes  the  real  part 

as  indicated.  For  all  practical  purposes,  this  is  the  case  for  Morlet  wavelets,  ip(t)  =  e”^'  e~^  ‘  which 
satisfy^  <  1/2  (cf,  [l]-[3],  Appendix  A),  Note  the  limits  of  the  integral  in  definition  (1 .1 ),  Some  references 
use  0  to  oo.  We  prefer  (1,1)  because  it  results  in  the  same  formulae  for  both  real  and  analytic  wavelets. 


where  the  are  tlie  so-called  duals  of  the  tp(2~'i  —  n).  Nevertheless,  under  propitious 
circumstajrces  the  wavelets  themselves  offer  a  good  approximation  to  the  duals.  One  way 
to  improve  this  approximation  is  to  add  voices;  i.e.,  discretize  a  at  smaller  intervals  (a  =  7', 
where  1  <  7  <  2)  [3].  This  technique  has  been  highly  successful  in  many  applications  [5]; 
however,  as  hinted  above,  it  is  not  sufficient  to  provide  a  good  inverse  for  narrowband 
wavelets  [6], 

If  we  are  to  study  the  situation,  we  must  be  able  to  define  the  duals  or  at  least  make 
their  roles  more  precise.  In  this  introduction,  to  .simplify  the  presentation,  we  avoid  defining 
the  functions  themselves.  Rather,  we  broach  the  subject  from  the  point  of  view  of  matrices 
and  their  pseudo-inverses.  There  is  an  additional  reason  for  doing  this.  Namely,  these 
concepts  are  particularly  appropriate  to  discrete,  finite  implementations.  Moreover,  as 
with  the  wavelets  themselves,  there  is  rarely  a  need  to  actually  compute  the  dual  functions; 
we  are  only  interested  in  the  wavelet  coefficients  and  the  inverse  algorithm.  Let  Wi^n  he 
the  decimated  discrete  wavelet  transform,  that  is,  the  discrete  analogue  of  1^(2*, 2'n).  It 
is  defined  by 


s.-t-i  -  i  if  *  s.)  (1.8a) 

w,  =  ^  *  s.  (1.8ft) 

where  [so]„  =  s(n)  is  the  discrete  signal,  /  and  g  are  discrete  filters,  and  J.  is  the  decimation 
(downsampling)  operator 


[ihU  =  h2n  (1-9) 

One  should  visualize  g,,  =  y^{~n)  as  the  sampled  wavelet  and  /  as  a,  somewhat  arbitrary, 
interpolation  filter  [2],  We  shall  not  dwell  on  regularity  questions  or  the  choice  of  the 
filters;  however,  we  remark  that  /  must  at  least  satisfy  the  condition 


yjn^V2.  (1.10) 


The  reader  is  referred  to  the  literature  [2]  -[3].  The  recursion  (1.8)  is  easily  implemented 
by  the  filter  bank  of  figure  1.  Ultimately  our  inverses  shall  be  based  on  the  undecimated 
discrete  wavelet  trairsform  (DWT).  It  corresponds  to  W(2',n)  and  is  illustrated  in  figure 
2a.  Denoted  wj,,  this  transform  satisfies  v/J.-,,  =  [2].  It  is  more  appropriate  to  most 

applications  of  nonorthogonal  wavelets  where,  as  we  have  mentioned,  output  resolution  is 
an  issue.  Moreover,  the  undecimated  transform  provides  a  better  basis  for  recovering  the 
signal  s(n). 


W2 


Figure  1.  A  wavelet  filter  bank  structure  implementing  the  decimated  discrete  wavelet  trans¬ 
form,  The  down-arrow  indicates  decimation.  The  output  of  the  transform  is  the  family  of  vectors 
Wi,  forming  the  two  parameter  transform  w,-,n  in  the  scale-time  plane. 


.1  D'  f  i 
s  - >-  s 


D'  g 


w 


(a)  DWT 


-  P'f  s' 


I 


D‘g 


w 


(b)  DVrr’ 


Figure  2.  One  stage  of  (a)  the  undecimated  discrete  wavelet  transform,  and  (b)  of  the  inverse 
discrete  wavelet  transform,  s®  is  the  input  signal  and  the  ’’initial”  stage  of  DWT~^  is  determined 
by  D7  is  the  filter  obtained  from  /  by  inserting  2'  —  1  zeros  between  its  elements. 


Suppose  there  are  M  stages^  to  the  DWT;  i.e.,  z  =  0, . . . ,  Af  —  1.  Then,  the  DWT 
maps  IR  into 

s  — *  {w  :  r  =  0, . , . ,  M  —  1  ;  s'"  )  .  (I  ll) 

Under  very  weEik  regularity  conditions,  the  norm  of  s'’^  goes  to  zero  as  M  goes  to  oo  (see 
Appendix  B).  However,  in  order  for  the  transformation  to  be  nonsingular  for  finite  M, 
one  must  include  the  smoothed  signal  as  well  as  the  wavelet  coefficients  w'.  This 
transformation  is  clearly  linear  (although  not  time  invariant  in  the  decimated  case)  and 
may  be  represented  by  a  matrix  A.  If  /  and  g  are  finite,  the  transformation  is  locally 
finite  (each  row  of  A  is  zero  except  for  a  section  of  fixed  length);  but,  the  matrix  itself 

^That  is,  so-called  scales  or  octaves.  In  general,  one  will  also  have  several  voices.  See  Section  IV 


is  infinite  because  the  convolution  acts  on  arbitrarily  long  signals.  Also,  the  image  of 
the  transformation  is  a  proper  subset  of  so  that  a  true  inverse  does  not  exist. 

Nonetheless,  it  is  generally  injective  so  that  A^A  is  nonsingular,  and  we  may  ’’invert” 
objects  in  the  range  space  ’  by  usixig  the  pseudo- inverse 

P  =  (AU)-’A^  (1.12) 

w’here  [A^],j  =  Aj^.  It  is  trivial  to  verify  that  PA  =  I.  On  the  other  hand,  AP  ^  I. 

There  are  many  other  matrices  Q  which  have  the  property  QA  =  I.  In  fact,  if  B  is 
anj'  matrix  such  that  BA  is  nonsingular,  then 

Q  =  (BA)-’B  (1.13) 


satisfies  .such  a  relation.  The  pseudo-inverse  is  unique  inasmuch  as  it  is  an  orthogonal 
projection  of  onto  the  image  of  the  DWT.  For  example,  suppose  one  were  to  filter 

the  wavelet  transform  (possibly  stepping  outside  of  its  range)  and  invert  it 
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ftiu 


h  €  1R^^+’ 


Q 


DWT 


Ar 


Then,  in  general,  Ar  ^  h,  but  the  pseudo-inverse  minimizes  the  squared  error  ||h  — Ar|JJ  = 
llh-AQh|p  among  all  possible  Q.  Note  that  this  holds  for  non-Euclidean  metrics  on  JR™'*'* 
provided  A^  is  replaced  by  the  adjoint  transformation  relative  to  the  metric.  Although 
the  pseudo-inverse  is  a  natural  choice,  other  issues  are  involved.  Often  the  error  criterion 
is  not  a  metric  on  IR^"''’,  or  computational  complexity  is  an  important  consideration. 

An  alternative  approach  is  to  invert  a  single  stage  of  the  DWT,  thereby  inverting  the 
entire  transform.  For  the  undecimated  DWT,  it  suffices  to  find  filters  /  and  g  such  that 


1  *  f  +  g*9  ~  ^  ,  (1-14) 

where  [5],n  =  <^o,fn  is  the  Kronecker  delta.  More  precisely,  if  D*  is  the  operator  which 
inserts  2’  —  1  zeros  between  the  elements  of  a  vector,  then  D'(/  +  <^)  =  D'/  +  D'^.  It  follows 
that  (1.14)  is  equivalent  to 


/Tkl  t\ 

V*-'  J  ) 


/TTll  C\  I  /I 

V*-'  1  )-r\i 


'  y; 


yij  g)  ^  o 


(1.15) 


Thus,  under  (1.14),  the  inverse  filterbank  of  figure  2b  is  a  left  inverse  for  the  DWT,  that 
is,  for  figure  2a. 

These  three  points  of  view  (dual  expansions,  left  or  pseudo-inverses,  and  filter  banks) 
begin  to  merge  if  one  considers  a  more  general  class  of  filters  than  those  satisfying  (1.14). 
The  usual  inversion  procedure  for  nonorthogonal  wavelets  is  to  approximate  (1.7)  with 


s(t)  Rt  ^  w,,„V^,,„(t) 

n,  i=0 . M  — 1 


constant 


•  i7e  ^  w, ,„!/>, >(<)  . 

n,  — 1 


I 


» 


$ 


(1.16) 


The  approximations  involved  in  the  traditional  frame  inverse  may  then  be  reinterpreted 
from  a  filter  bank  perspective.  The  finite  sum  (i.e,,  M  <  oo)  is  equivalent  to  ignoring  s^, 
while  the  use  of  the  wavelets  instf^d  of  their  duaJs  is  tantamount  to  setting  A^A  =  I  in 
(1.12).  Furthermore,  if  we  choose  /  and  g  to  be  /„  =  and  g,^  =  then  the  inver.se 
filter  bank  of  figure  2b  computes  the  adjoint  transformation  A^  (c.f.  Section  V).  Thus,  the 
inverse  filterbank  provides  a  unifying  framework  from  which  to  study  the  variou.s  inverses. 
In  fact,  we  shall  adopt  the  terminology  inverse  discrete  wavelet  transform  or  DWT~^  [G] 
to  refer  to  any  inverse  filter  bank  whether  or  not  the  filters  satisfy  (1.14).  To  be  useful, 
the  filters  should  be  chosen  to  provide  an  approximate  left  inverse  to  DWT;  however,  in 
general,  an  exact  left  inverse  will  only  be  achieved  with  iteration. 

There  are,  of  course,  tradeoffs  which  accompany  the  various  approximations.  Some 
of  the  negative  aspects  are  (a)  filters  satisfying  equation  (1.14)  are  usually  prohibitively 
long  (e.g.,  Morlet  wavelets),  (b)  these  filters  generally  provide  inverses  which  are  not  the 
pseudo-inverse,  and  (c)  although  AW  is  readily  computed,  (A^ A)"^AKv  is  not.  We  find, 
however,  in  the  context  of  inverse  filter  banks,  that  these  problems  may  be  successfully 
evaded.  First,  let  us  note  that  if  a  true  left  inverse  is  desired,  it  may  be  computed  at  a  eery 
moderate  computational  cost  simply  by  iterating  (so-called  Neumann  inverse)  the  forward 
and  reverse  transforms  of  figure  2.  In  fact,  computationally,  this  technique  seems  much 
more  effective  that  than  using  long  filters  to  achieve  comparable  accuracy. 

Secondly,  a  number  of  good  approximate  i’  verses  exist  which  are  sufficient  for  many 
applications  without  the  above  iteration.  For  example,  a  substantial  improvement  over 
(1.16)  is  obtained,  at  essentially  no  cost,  by  including  s^.  Also,  just  as  the  inclusion  of 
voices  provides  redundancy  creating  a  tighter  frame  and  therefore  a  better  inverse  [3],  the 
use  of  undc^cimated  wavelets  enables  us  to  find  inverses  which  perform  better  than  (1.16) 
with  considerably  less  computation.  ®  How  do  we  go  about  finding  these  inverses?  One 
very  effective  method  which  also  lends  considerable  insight  is  to  mimic  the  continuous 
case.  A  more  direct  approach  is  to  try  to  invert  a  single  stage;  that  is,  approximate  (1-14) 
using  filters  which  are  not  too  long  or  satisfy  suitable  constraints.  One  must  take  care, 
however,  to  insure  that  the  error  does  not  increase  unduly  as  the  inverse  passes  up  the 
filter  bank.  It  is  certainly  possible,  to  find  cases  in  which  some  degradation  of  the  single 
stage  approximation  actually  improves  the  performance  over  a  cascade  of  several  stages. 
It  may  also  be  preferable  to  use  a  coarse  approximation  and  then  iterate.  We  can  argue, 
however,  that  if  the  inverse  is  to  be  applied  to  an  unspecified  class  of  signals  and  operate 
with  little  or  no  iteration,  then,  a  single  iteration  of  a  one-octave  transform  must  give  a 
good  inverse;  i.e.,  the  single  stage  inverse  must  be  acceptable. 


“^In  general,  any  nonsingular  inverse  filterbank  is  a  linear  transformation  B  vvhich  may  be  iterated  in 
conjunction  with  the  forward  transform  to  compute  a  left  inverse  Q  defined  by  (1.13). 

^One  even  finds  that  the  duals  to  the  undecimated  wavelets  may  be  computed  (if  for  some  strange 
reason  one  wishes  to  do  so)  more  efficiently  since  they  are  translation  invariant  rather  than  scale  invariant 
as  in  the  decimated  case  (cf.  Section  III). 
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II.  DEFINITIONS  AND  NOTATION 


A  major  notatioiial  distinction  is  the  use  of  superscripts  for  undecimatod  entities 
versus  subscripts  when  they  are  decimated.  We  define  the  standard  families  of  wavelet 
functions  generated  from  the  mother  wavelet  i/j  by 


*  (^) 

’Pi.nit)  =  ■ 


The  corresponding  sampled  continuous  wavelet  transforms  are  given  by 


=  J  rp\(T)  dr  =  <s,rpl> 

=  j  s{t)  dr  =  <s,V^.,„>  . 


(2,1a) 


(2.U) 


(2.2o) 


{2.2b) 


This  notation  differs  slightly  from  [3]  which  uses  V’"'”;  however,  it  enables  us  to  treat  the 
discrete  transforms  as  vectors;  i.e.,  wj,  =  [w*]„.  Unless  specified  to  be  otherwise,  norms 
are  assumed  to  be  for  functions,  l|t/i||^  =  and  for  vectors,  ||w||^  —  X^„(wn)^- 

Finally,  let  us  remark  that  additional  voices  are  indicated  by  and  V’tt.n,  which  are 
equivalent  to  (2.1a)  and  (2.1b)  respectively  with  2'  replaced  by  2*  •  2'’^^  where  L  is  the 
number  of  voices.  More  formal  definitions  are.  postponed  to  Section  IV  where  voices  will 
be  treated  in  detail. 

For  our  purposes,  Wi  and  W‘  serve  mainly  as  an  aids  in  the  interpretation  of  the 
corresponding  discrete  wavelet  transforms.  These  discrete  transforms,  w;,,,  and  wj,,  may 
be  computed  according  to  the  diagrams  of  figures  1  and  2a  respectively.  The  former  is 
defined  by  equations  (1.8);  the  latter  is  defined  eis  follows;  Let  D‘  be  the  dilation  operator 
that  inserts  2'  —  1  zeros  between  each  of  the  elements  of  a  vector;  i.e.. 


ID'/kwr  = 


for  r  =  0 
otherwise 


For  a  given  lowpass  filter  /  and  highpass  filter  g,  the  undecimated  DWT  is  determiixed  by 
the  following  recursion  [2] 

g.+i  ^  d7*s'  (2.4a) 

w*  =  D7*s'  (2.46) 

with  [s®J„  =  s(n),  the  discretized  signal.  The  filter  g  is  generally  determined  by  sampling 
the  mother  wavelet;  that  is, 

[?]n  =  w{~^)  (2-5) 


I 


The  relationship  between  the  two  transforms  may  be  expressed  by 


Wi.„  =  W,i, 


We  make  the  standard  observation  that  the  are  translation  invariant 


V’n  =  [V’‘]n  =  >  (2.7a) 


while  the  V’i.n  not 


V«4V  y'tjrt  UO. 

V’l.n  =  ^  •  (2.7fc) 

Similarly,  let  be  the  operation  of  translation  by  m;  i.e., 

[TmS]n  =  S„_m  ,  (2.8) 

and  introduce  a  vector  argument  s  to  w'  to  indicate  dependence  on  the  signal.  It  follows 
immediately  from  eq\.iations  (2.4)  that  the  undecimated  wavelet  transform  is  translation 
invariant,  i.e., 

w|.+m(s)  =  w;,(T_„,s)  .  (2.9) 

On  the  other  hand,  as  a  consequence  of  the  downsampling  in  (1.8),  the  w,  „  are  not. 

We  shall  use  a  tilde  to  denote  the  dual  wavelets.  This  operator 


_  (^1,.  \~ 
<pt,n  —  \Yi,n) 


(2.10) 


depends  on  the  structure  of  the  entire  family  of  {V’i.m},  not  just  on  the  single  function 

As  a  consequence  the  translates  (dilates)  of  the  dual  wavelets  are  not  necessarily  the 

duals  of  translations  (dilations)  of  the  wavelets.  That  is,  the  dual  operator  does  not  always 
commute  with  translation  and/or  scaling.  In  particular,  although  )/’o,o(t)  =  ~  V’(0> 

we  have  the  bizarre  effect  that  _ 

V'o,(i(0  7^  V-SIO  1  (2-11) 

and,  in  fact,  makes  no  sense  at  all.  We  stretch  the  notation  by  also  using  tilde  to 
denote  the  niters  of  the  inverse  niter  bank,  e.g.,  f  »,nd  y.  Tlie  idea  is  that,  just  as  the 
forward  filters  are  associated  with  the  wavelets,  the  inverse  filters  arc,  at  least  conceptually, 
associated  with  the  dual  wavelets. 

The  operator  f  has  several  (consistent)  interpretations  depending  on  the  context.  For 
matrices,  it  denotes  the  adjoint 

[A%  =  [Xy  .  (2.12a) 

For  vectors,  it  is  the  complex  time  reversal 


[/'].  =  [/]-.  , 


(2.126) 
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which  is  consistent  with  the  adjoint  of  the  matrix  corresponding  to  convolution  by  /  ; 

that  is,  with  respect  to  a  matrix  [F]i,j  =  /,_j  where  f  *  s  =  Fs.  More  generally, 
represents  the  adjoint  of  an  operator  in  a  Hilbert  space  with  inner  product  <• ,  •> 

<A^x  ,x>  =  <x,  Ax>,  (2.12c) 

which,  under  the  Euclidean  metric,  agrees  with  the  matrix  adjoint.  Finally,  we  shall  have 
some  occasion  to  use  the  z-transform.  It  is  defined,  on  the  unit  circle,  by 

/»  =  H/ne''"”.  (2-13) 

n 

which  agrees  with  the  Fourier  transform  defined  in  (1.2).  Note  that 

(D7).(w)  =  /,(2-u,)  .  (2.14) 

» 


I 


I 


i 


III.  PROPERTIES  OP  THE  DUALS 


In  this  section,  we  define  the  dual  wavelets,  relate  them  to  the  pseudo-inverse,  and  de¬ 
rive  their  properties  under  translation  and  scaling.  We  demonstrate  that,  unlike  the  duals 
of  the  decimated  wavelets  (cf.  [3]),  undeciinated  wavelets  have  duals  which  are  translation 
invariant.  The  unmathematical  reader  may,  without  loss  of  continuity,  shun  all  but  subsec¬ 
tion  III.B  and  simply  consider  the  dual  wavelets  to  be  an  arbitrary  family  of  functions 

that  provide  the  inverse  expansion  s{t)  =  ^  Wi  ,n^,,nW  =  V’.,n(0-  With 

a  little  more  effort,  one  might  wish  to  keep  in  mind  the  following  discrete  interpretation  (cf. 
Appendix  C):  Let  A  be  a  matrix  representing  the  di.screte  wavelet  transformation,  and  A^ 
its  adjoint.  Then,  the  rows  of  A  correspond  to  the  frame  vectors  with  component 
[V'i,n]A:  =  and  thc  dual  vectors  0,_n  correspond  to  the  columns  of  (A'^A)“^A^ 

(cf.  the  pseudo-inverse  of  equation  (1.12)).  Note  that  A  depends  on  the  entire  family  of 
wavelets  so  that,  V-'i.n  has  a  nonlinear  dependence  on  rpi^n-  The  orthonormal  case,  where 
A^A  is  the  identity  so  that  i/h,n  --  is  an  exception. 

For  continuous  signals  and  real  wavelets,  the  dual  wavelets  are  given  by® 

0i,n  =  V’]t,n  i  (3.1) 

where  the  wavelet  transform  is  defined  by  (Ajec =  <s,  The  analogous 

expression  is  u.sed  to  define  thc  dua,ls  j/jJ,.  Different  metrics  on  the  image  space  produce  dif¬ 
ferent  adjoints  and  hence  different  families  of  duals  and  different  left  inverses.  In  particular, 
the  standard  frame  formulation  uses  a  re.solution  of  the  identity 

•S  —  <^3,  (3-2) 

i,n 


based  on  thc  Euclidean  metric  for  decimated  wavelets.  That  is,  it  is  the  pseudo- inverse 
(1,12)  where  is  defined  relative  to  the  metric  l|W(jetll'^  =■•  and  4>i,n  by 

(3.1).  This  metric  corresponds  to  energy  since  is  power  per  Hz  and  the  cells  at 

octave  i  have  area  (2'  secj(l/2'  Hz)  —  1  (see  Section  IV).  Similarly,  the  ciiergy  metric  fur 
undeciinated  wavelets  is  ||w“”‘^®‘^|l*  =  „  2'*  IwJJ"'',  so  that,  for  example,  if  thc  decimated 

wavelets  arc  orthonormal,  the  undecimated  dual  wavelets  satisfy  y')]^  = 
contrast  to  =  V-'en-  The  energy  metric  is  intuitively  ])leasing,  and  for  decimated 
wavelets  its  pseudo-inverse  is  an  orthogonal  projection,  so  this  left  inverse  is  often  the  one 
of  choice.  However,  one  should  not  arbitrarily  dismiss  alternative  metrics.  First,  as  we  shall 

®]n  the  case  of  complex  analytic  wavelets,  a  complete  set  for  1/*(1R)  i.s  given  by  {V-'ji  V^)}  is, 

strictly  speaking,  the  frairie  expansion  involves  two  niothcr  wavelets  rjj  and  yp.  Similarly,  the  family  of 
duals  is  {ypj,  V' j  }  ['ll-  Alternatively,  for  real  signals,  one  may  take  thc  real  j>art  as  in  (1.7).  This  aspect  is 
ignored  in  this  section  in  order  to  sirriplil'y  the  presentation. 


see  ill  succeeding  sections,  other  left  inverses  have  considerable  coinputational/numerical 
advantages.  Furthermore,  it  is  not  a  priori  clear  that  alternative  metrics  may  not  enhance 
one’s  signal  processing,  at  least  in  particular  situations.  In  fact,  every  inverse  of  the  form 
(1.13)  corresponds  to  a  pseudo-inverse  with  respect  to  an  appropriate  metric. 

We  make  this  more  precise.  The  pseudo-inverse  may  be  characterized  as  is  the  unique 
left  inverse  such  that  its  null  space  is  orthogonal  to  the  range  of  the  transform.  We  note 
that,  provided  BA  is  nonsingular,  Q  =  (BA)~*B  is  a  left  inverse  of  A.  Also,  since  BA 
is  nonsingular,  the  null  space  of  Q  coincides  with  the  null  space  of  B,  and  it  is  easy  to  see 
that  the  reinge  of  A  and  the  null  space  of  B  are  disjoint  and  span  the  embedding  space 
(i.e.,  time-scale  space).  Thus,  any  vector  w  in  this  space  may  be  written  uniquely  as  the 
sum  of  vectors  in  these  two  subspaces.  For  Q  to  be  the  pseudo-inverse,  it  suffices  to  chose  a 
metric  for  which  these  two  spaces  are  orthogonal.  In  fact,  such  a  metric  is  given  explicitly 

by 

<w,  w>B  =  1|A(BA)-‘Bw1|2  d-  11(1- A(BA)'’B)wf  .  (3.3) 

A.  Invariance 

In  [3]  it  is  shown  that  the  duals  of  the  decimated  wavelets  are  scale  invariant  in  the 
sense  that  V’«,n  =  but  that  they  are  not  invariant  under  translation,  i.e., 

Vd,n  V’i,o(^  ~  ^■)-  In  this  section,  we  show  that  for  undeciiiiated  wavelets  the  reverse 
is  true.  That  is,  their  duals  are  translation  invariant  but  not  scale  invariant.  This  is  of 
theoretical  interest  and  also  has  a  certain  degree  of  practical  impact.  It  implies  that  the 
undecimated  duals  may  be  all  computed  by  determining  one  function  for  each  octave  and 
then  translating.  Since  there  are  generally  many  more  time  points  than  there  are  octaves, 
this,  property  appear.s  to  be  more  useful  than  scale  invariance. 

The  for  ward  undecimated  waveiet  transformation  A  may  be  represented  by 

M.s];,  =  <s,  rn>  ■  (3.4) 

Wc  then  have  (with  the  energy  metric) 

<A^As,  r>  =  <As,  Ar> 

-  <.9,  !/';,>  ,  (3.5) 

i,n 

which  implies  that 

am  =  532-‘<.,  Vh\>V.’^  .  (3.6) 

1,71 

Next,  let  us  define  the  scale/tirne  translation  operator 
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# 


Its  adjoint  is  given  by 


(3.8) 


which  follows  from 


J s(t))  r(t)  dt  =  J  r(t)  dt  =  J  s{t)  r^ra  +  b)  y/a  dr  .  (3.9) 


Note,  that  we  also  have 


Tr(“*)  Tr(‘^)  _ 

~  ^b'+a'b  • 


(3.10) 


We  are  particularly  concerned  with  the  family 

Ui,  ^  up  .  (3.11) 

If  A  commutes  with  then  the  family  of  duals  is  translation  invariant;  if  it  commutes 
with  Uq,  then  the  family  is  scale  invariant.  For  example,  commutation  of  A"^  A  with 
implies  commutation  of  (.A^  .4)~^  with  and 

^(<)  -  (-4U)-Vj,  =  {A^AVulr,  =  UliA^Ap^l  -  ^•(t-n)  .  (3.12) 

Proposition  S.l:  Let  ip  be  any  admissible  wavelet.  Then, 

(i)  Under  any  transla.tion  invariant  metric,  the  family  of  dual  wavelets  for  the  undeci¬ 
mated  discrete  wavelet  transform  is  time  invariant  but  not  scale  invariant;  i.e., 

^^(t)  =  ^‘(t-n);  ^•(t)  ^  2-‘/^^o(t/2‘)  .  (3.13) 

(ii)  Under  any  scale  invariant  metric,  the  family  of  dual  wavelets  for  the  decimated  discrete 
wavelet  transform  is  scale  invariant  but  not  time  invariant;  i.e., 

pijt)  =  2-/^:^,(f/2*)  (3.14) 

We  remark  that,  in  l°°,  a  translation  invariant  metric  <x,  y>g  =  Ylijnm  ^h9nmyL  one 
for  which  Scale  invariance  is  defined  analogously.  For  notational  clarity, 

we  only  prove  (i)  under  the  energy  metric. 

Proof:  From  (3.6),  we  compute 

U>,  A<A  =  £2-  <-.U^4,>  C/'.C' V.  =  !;<■,  ,  (3.18) 


A i  A  rri 
U  - 


\  ’’  rt  —  i  _  TTltm'  I.  Trt  » 

2^  ^  .  I'm  '.'n  VX.'.i  V 


>,n 

=  vii, 

i,n 

=  <  ■.  c;'-, 


(3.16) 


Even  if  we  accept  translation  invariance  up  to  a  scalar  factor  2“'^,  to  equate  (3.15)  and 
(3.14)  we  still  need  to  make  the  substitution,  n'  =  2~^(n  —  m).  This  is  only  possible  for 
j  =  0  since  n'  must  be  an  integer.  In  fact,  for  general  if),  equality  occurs  only  for  and 
not  for  Ul .  On  the  other  hand,  for  decimated  wavelets,  the  operator  A^^^Adtc  takes  the 
form  „  <  •,  it  is  easy  to  confirm  that  this  operator  commutes 

with  1/q  ,  but  not  with  U^. 

B.  The  Neumann  Inverse 

Often,  given  a  discrete  wavelet  transform,  A  and  an  approximate  inverse  B  (for  ex¬ 
ample,  an  inverse  filter  bank  DWT~^),  we  would  like  to  compute  an  exact  left  inverse,  in 
particular  that  of  equation  (1.13).  This  may  be  done  by  using  the  Neumann  inverse  for 
linear  operators^ 

OO 

5-‘  =  ,<  (3.17) 

Jt=0 

where  S  is  a  Hilbert  space  operator  and  p  is  sufficiently  small  to  insure  convergence;  that 
is,  (||(i  —  /i5)ll  <  1  [3],  [7].  Equation  (1-13)  may  then  be  written 

OO 

(BAr^B  =  ^  Y,^I-(jiBAfB  .  (3.18) 

fe=0 

The  summands  in  (3.18)  can  be  computed  by  the  recursion  ejt  =  e*_i  —  fiB  Aei-\.  Here, 
A  and  B  are  the  transforms  DWT  and  DWT~^  ^  implemented  by  filter  banks.  One  may 
use  the  following  algorithm  to  cora,pute  a  true  left  inverse  as  accurately  as  desired. 

Algorithm  S.l 

(0)  xo  =  {DWT)-^w 
Co  =  Xo 

(1)  €k  =  tk-i  -  fi  DWT^^  o  DWT  u-i 

(2)  xjt  =  xjt_i  +  ek 

(3)  leftinverse{k)  —  fiXk 
go  to  (1) 

A  computational  drawback  that  becomes  apparent,  up'-n  applying  Algorithm  3.1,  is 
the  huge  increase  in  the  signal’s  support  (length)  under  the  mapping  DWT~^  o  DWT. 
In  fact,  if  a  filter  g  has  length  r  then  D‘j  has  length  (2'  —  l)(r  —  1)  -t-  r,  and  (D'g)  *  s* 
adds  (2’  —  l)(r  —  1)  points  to  the  support  of  s*.  There  is  an  additional  effect  due  to  D'/, 
although  for  Morlet  wavelets  g  is  much  longer  than  /  and,  hence,  dominates  tl.\e  problem. 

^This  formula  is  eiisily  understood,  if  we  replace  the  linear  operator  by  a  scalar  s,  let  x  =  1  —  pa, 
and  consider  the  expansion  =  p/(l  — x)  =  p 


The  same  comments  hold  true  for  DWT~^.  For  example,  for  five  octaves  and  r  =  50,  the 
time  series  DWT  o  DWT~^w  will  be  about  3000  samples  (2  x  31  x  49  )  longer  than  w. 
Each  iteration  of  Algorithm  3.1  adds  such  an  amount  to  the  support  of  xjt-  Clearly,  this 
is  unsupportable.  Fortunately,  the  obvious  heuristic  solution,  clipping  the  support  of  the 
iterations,  works  and  can  be  justified.  Let  be  the  operator 


(^w  y)’„ 


j  yj,,  for  (i,n)  in  the  support  of  w 
\  0,  otherv/ise 


(3.19) 


We  simply  replace  DWT~^  in  step  (1)  of  Algorithm  3.1  by  Tiw  o  DWT~^ .  That  is,  we 
set  DWT  €k  to  zero  outside  of  the  support  of  w.  This  is  equivalent  to  replacing  B  in  (3.18) 
by  B'  =  Hy,  B,  and  leads  to  a  valid  left  inverse  provided  p  is  small  enough  and  B'  A  is 
nonsingular.  Because  of  the  w-dependent  clipping,  this  inverse  is  actually  nonlinear.  Of 
course,  one  could  always  create  a  similar,  but  linear,  operator  by  restricting  oneself  to  a 
family  of  signals  with  fixed  support, 


IV.  APPROXIMATIONS  TO  CONTINUOUS  INVERSES 


Of  course,  (1.5)  is  not  the  only  continuous  inverse.  As  with  the  wavelet  scries,  there 
is  an  infinity  of  left  inverses  for  the  continuous  wavelet  transform.  In  this  section,  we 
examine  two  of  these.  The  results  shall  be  used  in  Section  V  as  guidelines  for  inversion 
where  they  provide  discrete  metrics  and  appropriate  filter  normalizations  for  the  DWT~^ . 
In  the  presence  of  voices  these  normalizations  become  nontrivial  and  extremely  confusing; 
however,  they  are  essential  for  correct  implementation  of  the  inverse  transform. 


A.  D oublt-initgral  inverse  (resolution  of  identity) 

We  begin  with  (1.5),  which  for  convenience  we  rewrite  as 

(41) 

where 

Recall  that  s  is  real,  rj)  must  be  real  or  analytic  (V’(‘^)  =  0  for  w  <  0),  a  is  nonnegative,  and 
Re  denotes  the  real  part.  The  wavelet  transforms  themselves  are  generally  complex.  For 
simplicity,  we  shall  sometimes  drop  the  Re.  Also,  observe  that  a  factor  l/\/a  is  included  in 
the  definition  of  the  wavelet  (0-  a  result  the  norm  of  is  independent  of  a  and 
b.  It  often  helps  conceptually  to  normalize  the  mother  wavelet,  =  1.  Alternatively, 
the  norm  of  t/'  may  be  absorbed  in  the  constant  (cf.  equation  (1.1)). 

Let  us  discretize  the  integral  and  compare  it  to  the  sum  (1.7).  The  coefficients  are 
given  for  for  all  integral  values  of  b  in  the  undecimated  case  and  for  6  =  2'  in  the  decimated 
case.  This  results  in  d6  «  A6  =  1  and  dh  r?  2'  respectively.  When  the  scales  are  output  by 
octaves,  a  takes  on  the  values  a  —  2'  for  integer  i.  Thus,  da/ a  =  d(lna)  Rs  A(zln2)  =  In  2, 
Substituting  these  expressions  into  (4.1)  yields 

s(i)  =5  (4.3a) 

^  i,n 

for  the  undecimated  wavelets,  and 

s{t)  sj  (4.36) 

L'lJ, 

^  »,n 


for  the  decimated  v/avelets  since  A6  =  2'  cancels  the  factor  in  the  denominator.  Comparing 
these  with  (1.7)  and  its  undecimated  equivalent,  we  see  that 


V4 


21n2  1 
2' 


(4.4a) 


(4.46) 


2 In 2  /.  .1^, 

^i,n  «  -TT—  V>i,n  •  (4-46) 

Note  that  (4.4b)  agree.s  with  (1.16)  and  that  both  the  expressions  (4.4a)  and  (4.4b)  are 
exact  when  the  form  an  orthonormal  basis  [3]. 

However,  we  are  interested  in  iionorthogonal  wavelets,  which,  in  general,  require  a 
finer  output  grid,  that  is,  voices.  If  each  octave  is  divided  into  L  voices,  we  have 

a  =  2'y'’  =  ,  (4.5) 

where  i  =  0, . . .  M  —  1;  u  =  0, . . .  L  —  1;  and 


7  ^  2^/^ 


Substituting  (4.5)  into  (4.2)  yields 


^/¥r ( 2-r  )  ’ 


with  a  similar  expression,  =  2'‘/^7~*'/'‘' V’(2“*7“*’t  —  n),  for  decimated  wavelets. 

Voices  for  the  discrete  transform  are  usually  computed  by  replicating  the  filter  bank  with 
a  new  mother  wavelet 


‘■(0 


for  each  voice.  Althorigh  in  the  undecimated  case  this  remains  equivalent  to  (4.7),  in  the 
decimated  case,  it  implements  a  DWT  corresponding  to  the  hybrid  family  of  wavelets 


Decimated  voices  are  rarely  implemented,  but  when  they  are,  one  invariably  uses 
rather  than  V'jTJ.n- 

Proceeding  on  the  basis  of  (4.7)  and  (4.9),  we  retain  the  output  grids  A6  =  1  for 
undecimated  wavelets  and  A6  =  2'  for  decimated  wavelets,  while  incrementing  a  in  (4.5) 

ViTf  1  ^7 4"  V«  4  /  jVk  »»  n  r\'^  n  •  T  «  r*  « •Tro  1  4  w\  ^  rv  •  i  ^  —  1  *  i  «  n 

•-/Jr  V  .At  U*  V/X  K/  »jj  U  *«J  \.y\^  U  A  V  U  1.1V/  tXiXi  V/V'VU'VVy^  X  •  Vy .  ^  VV/  “  i.  «  J.4.XVt.O) 

d(lna)  w  A(z  +  u/X)ln2  =  \n2/L.  Substituting  these  expressions  into  (4.1)  yields 


s{t)  ^  V"  VV;"  ~V’;'’(<)  • 


(4.10a) 


for  the  undecimated  wavelets,  and 


(4.106) 
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(4.116) 


for  decimated  wavelets.  The  corresponding  duals  are 

^ 

It  should  be  empheisized  that  these  expressions  for  the  duals  are  mainly  intended  to  serve 
zis  conceptued  guidelines.  For  narrowband  wavelets,  the  quality  of  the  approximation  is 
certainly  in  question. 

The  normalizations  (4.7)-(4.9)  have  several  advantages.  The  squares  of  all  three  types 
of  wavelet  coefficients  represent  power  per  Hz.  This  follows  from  the  “Plancherel”  formula 
(cf.  Appendix  D) 

Energy  =  j\s{t)\'‘dt  =  A  jf“  ^ 

r  t  *1  n  ' 


(4.12) 


Equation  (4.7),  implies  that  the  bandwidth  of  is  proportional  to  l/(2'7'')  so  that  the 
measure  in  the  integral  (4.12)  is  in  units  of  seconds  ■  Hz.  This,  in  turn,  im-plies  that  the 
essentially  power/ Hz  [2].  The  same  remains  true  of  Wiv,n  and  since  they  are 
simply  nonuniform  samples  of  the  same  density  function;  that  is,  of  VV'"''(n).  Furthermore, 
this  normalization  may  be  simply  characterized  as  that  which  maintains  the  norms  of  the 
analyzing  wavelets: 

=  Uiv,u\\  =  =  ll’^ll  .  (4.13) 

A  different  density,  perhaps  more  in  the  spirit  of  scaling  functions,  i.s  power/scale 
where  the  units  of  scale  are  ln(a).  This  normalization,  used  for  example  in  [8],  results  in 
analyzing  wavelets  of  the  form 


■Z*7"’  \  ■Z‘7''  y  a  \  a  J 


(4.14) 


Another  normalization,  the  symmetric  one  of  [3]  approximates  self-dual  wavelets;  i.e., 
V’/n  ~  const  •  V’,“,n  where  the  constant  is  independent  of  i  and  v.  This  is  achieved  by 
including  the  square  root  of  the  coefficient  appearing  in  (4.11b)  in  the  normalization  of  i/>". 
Thus,  ^  gets  replaced  by  is,  by  2''*/^7~''^(2"'7“"t  —  7“"n).  However, 

this  symmetric  normalization,  results  in  transforms  which  are  not  densities.  Their  squares 
are  snergy/cell  where  the  cell  size  varies  with  scale  and,  hence,  is  difficult  to  interpret. 
As  far  as  the  composition  DWT~^  o  DWT  is  concerned,  the  normalizing  factor  may  be 
arbitrarily  split  between  the  wavelets  and  their  duals.  However,  (a)  different  choices  will 


»  • 


17 


clijmgc  the  appearaiice/interpretation  of  the  wavelet  transform  image,  (b)  each  choice 
corresponds  to  a  different  left  inverse,®  and  (c)  the  implementations/approximations  of  the 
various  in'^'tses  may  exhibit  significantly  different  behavior. 


B.  Single-integral  inverse  (Morlet  formula) 

A  second,  less  frequently  used  formula,  is  the  so-called  Morlet  inversion  formula  ( [8] , 
Appendix  D) 

where  either  if}  is  analytic  or  t/>(cu)  is  real, 


C 


\g> 


A 


r  ^ 

J-oo  1^1 


duj  , 


(4.16) 


and  Wf  =  defined  in  equation  (1.3).  For  Morlet  wavelets  ip{uj)  is  real,  and,  for 

practical  purposes,  ip(t)  is  analytic.  The  subscript  of  C\xi)  is  intended  to  reflect  the 
dependence  on  ip  in  (4.16)  as  opposed  to  its  modulus  squared  which  appears  in  (1.1).  The 
computational  advantage  of  a  single  integration  is  enticing,  and  we  shall  find  that  the 
corresponding  discrete  formula  is  quite  effective  when  applied  to  undecimatcd  wavelets. 

Discretizing,  as  in  the  previous  subsection,  we  obtain 


s{t) 


Re 


2  ln2 
Ci^  L 


V- 


-  wr 


(4.17) 


The  continuous  time  dependence  of  the  coefficients  IF,**'  in  (4.17)  prevents  a  direct  inter¬ 
pretation  as  a  dual  expansion;  however,  if  we  discretize  the  output  emd  define  the  vector 
6^  by  =  ^nnt)  then  we  may  rewrite  that  expression  as  s  w  where 


2  ln2  r 

tt>  — _ . 


(4.18) 


The  formulae  for  the  decimated  wavelets  are  essentially  the  same  since  there  is  no  integra¬ 
tion  over  b.  However,  when  we  discretize  the  time  by  setting  w^„^„  = 

»n  rM'#»QAnf  (\  a  nnf.  all  nrf.avpc  Thp 

^  - - - - -  —  - - - J  —  - - - - - — .  - 

lack  of  these  values  can  have  a  crippling  effect  on  the  inversion  approximation  (cf.  Section 

V). 

We  remark  that  both  of  the  inversion  formulae  that  we  have  presented  are  special 
cases  of  a  more  general  family  of  inverses  (cf.  Appendix  D) 


1  dn  r°°  _ 

»(*)  =  7=r- /  <».(!(''');>  w'')nc 

^  AB  J—oo  ®  J  ~oo 

®That  is,  the  duals  are  computed  under  different  metrics,  cf.  Section  III 


(4.19) 


where 


J —oo 


(4.20) 


Equation  (4.1)  is  derived  bj'  setting  ?/'^(t)  =  V'^(^)  =  V’(0t  while  (4.15)  follows  from 
=  ^(t).  Additional  manipulations  and  restrictions  on  ip  allow  one  to  replace  the 
integral  over  negative  a  by  twice  the  real  part. 


C.  Rtlc.tion  to  DWT 

Let  us  now  relate  some  of  the  above  approximations  to  the  discrete  wavelet  transform. 
We  examine  the  undccimated  DWT  of  figure  2a  invoking  the  approximation  wf,  It 

is  clear  that,  if  there  are  M  stages,  then  s*'  replaces  the  coefficients  W*  for  i  >  M.  If 
we  consider  s"  as  a  discretized  signal  and  take  its  wavelet  transform,  (4.12)  with  a  single 
voice  implies 

--■(rinp  .  (4-21) 


where  ||w’|p  =  l|wnl|^.  Substituting(4.21)  back  into  (4.12)  then  yields 


\mf 


21n2  ^  1 


E  + 


(4.22) 


When  its  approximations  hold,  this  equation  implies  a  relationship  between  the  normaliza¬ 
tions  of  the  filters  /  and  g.  Consider  the  output  of  a  one-stage  (M  =  1)  filterbank  where 


i.1.  _  j: _ i..  _ t  ; _ i-.  -  —  f  \t(t^  II„I|  _  1  »,r0  _  ..  _  f 

hue  uiaeicte  aigiidi  iS  an  iiiipuiac  luneutuii  5  —  %/•  yvc  *iavt.  |js||  . —  ”  —  J#)  aiiu  o  —  j 


giving 


^IIjIP  +  ^II/IP  «  1  ■ 


(4.23) 


Multiple  voices  are  easily  implemented  by  making  L  copies  of  the  filter  bank,  each 
with  a  different  filter  jr” 

(r )  ' 

Note  that  a.s  long  as  (4.24)  is  a  good  piecewise  approximation  to  ip,  we  will  also  have 
iiy^ii  ^  iiV'ii)  which  we  shall  always  assume  to  be  one.  The  general  expression  (4.12)  for 
the  total  energy,  with  voices,  takes  the  form 


v=0 


(4.25) 


while  (4.23)  becomes 


TEills”!'  +  5II/II'  “  1  ■ 
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(4.26) 


'■jL 


t 


It  is  not  difficult  to  show  that  (Appendix  A) 


L~\ 


1  1 
2-^'^  ~  o  n  -  9- 


t)=0 


2(1 


21n2 


(4.27) 


so  that  (4.26)  becomes 


J_  1 
^  2 


1  . 


(4.28) 


To  be  energy  preserving,  the  wavelets  must  not  be  too  redundant  or  too  sparse.  If  they 
exactly  fill  the  spectrum;  that  is,  if  const  ■  7~*’|<7”(a’)P  +  ^l/z(a^)|^  =  I5  then  the 

approximations  rire  exact  [9]. 


D.  Examples 

We  examine  the  quality  of  the  above  approximations  relative  to  our  prototype  non- 
orthogonal  wavelet  transform,  namely,  the  DWT  with  g  based  on  a  Morlet  wavelet  and  / 
an  a  trous  filter  (i.e.,  /2„  =  (5„o / >/2,  [2]).  The  Morlet  wavelet  of  norm  one  takes  the  form 
V’jv(t)  =  Its  v^'‘  voice,  7~*’/^V’m(7~"0(  has  a  relative  bandwidth  of 

P/'y'’-  Requiring  that  these  voices  cover  the  upper  half  of  the  spectrum  results  in  the 
following  condition,  which  relates  L  and  /3  (cf.  Appendix  A) 


L-i 


t;=0  ' 


TT 


(4.29a) 


or,  from  (4.27) 


L  > 


In  2 

ln(l-0.9jS) 


(4.296) 


We  remark  that  the  minimum  number  of  voices  satisfying  (4.29)  is  approximated,  conser¬ 
vatively,  by  T  =  1/13.  is,  of  course,  a  function  of  13.  It  takes  on  the  values  2.6,  2.26,  and 
2.12  for  /3  equal  to  0.5,  0.25,  and  0.125  respectively.  As  (3  approaches  0,  approaches 
2.0.  That  is,  the  set  g'’  acts  as  a  perfect  highpass  filter  (cf.  ((4.28)).  For  L=l,  equation 
(4.29b)  is  an  equality  for  (3  =  0.56.  That  value  of  implies  that  =  2.76  which  results 
in  21n2/C'^  —  .502.  For  L=  8,  we  get  [3  =  0.092,  =  2.09,  and,  using  ||j"||  =  1  along 

with  (4.27),  we  find  that  the  first  term  in  (4.26)  has  the  value  0.499. 

Two  typical  interpolating  filters  f  are  the  Lagrange  a  trous  filters  (1/2,  1,  1/2)  /  \f2 
and  (—1/16,  0,  9/16,  1,  9/16,  0,  —1/16)  /  \/2  which  have  squared  norms  of  0.75  and  0.82 
respectively.  Ideally,  we  would  like  the  lowpass  filter  /  and  the  set  of  highpass  filters  g'’  to 
split  the  energy  spectrum,  into  two  equal  parts  [2].  In  reality,  /  and  g  are  not  orthogonal 
so  that  some  energy  is  bound  to  leak  across  scales.  In  fact,  one-half  of  the  above  values 
are  0.375  and  0.41  respectively,  so  that  (l/2)||/|p  equals  0.5  only  very  approximately. 
This  can  be  interpreted  as  a  need  to  compensate  in  the  energy  expression  for  the  lack 


of  orthogonality  between  the  lowpass  and  highpass  filters.  It  would  seem  reasonable  to 
correct  this  discrepancy  by  actually  choosing  in  (4.26)  to  satisfy 


2  In  2  1 

~C^L  ^ 

^  V 


7" 


(4.30) 


where  ||y''||  =  1. 

As  another  interesting  example,  we  examine  the  Haar  wavelets.  Appropriate  DWT 
filters  are  /  =  (1,  l)/v/2  and  y  =  (— 1,  l)/\/2,  resulting  in  exact  wavelet  coefficients  (if  s° 
is  in  the  appropriate  subspace)  for  a  continuous  transformation  with  the  mother  wavelet 


r-l/v/2,  for-l<t<0 

V'3(i)  =  <  1/V^,  forO<t<l  •  (4-31) 

0  otherwise 


We  have  added  the  subscript  g  because  the  Haar  wavelet  is  more  commonly  taken  to  be 
\/2.  The  difference  springs  from  maintaining  consistency  with  the  relationship  gn  = 
7/>3(— n)  of  the  generalized  DWT,  and  is  reflected  by  the  fact  that  the  usual  orthonormal 
filter  bank’s  first  stage  output  is  interpreted  a.s  scale  1  rather  than  scale  0  (see  [2]).  A 
nontrivial  computation  yields  =  4  In  2.  Consequently,  the  first  coefficient  in  (4.23)  is 
(21n2)/(41n2)  =  1/2.  Since  the  norms  of  f  and  g  are  1,  (4.23)  holds  exactly.  This  exact 
result  is  somewhat  surprising  inasmuch  eis  was  computed  from  the  continuous  wavelet 
and  the  sums  in  the  formulae  were  obtained  by  approximating  an  integral. 


V.  THE  DISCRETE  INVERSE 


A  major  goal  of  this  section  is  to  develop  and  redefine  previous  concepts,  such  as  energy 
and  the  inverse  wavelet  transform,  entirely  in  discrete  terms.  We  shall  find  this  tantamount 
to  a  specification  of  filterbank  implementations  and  a  description  of  their  properties.  In 
keeping  with  this  spirit,  we  reserve  the  term  discrete  inverse  wavelet  transform  DWT~^ 
to  mean  precisely  the  family  of  transformations  computed  by  the  filter  bank  of  figure 
2b.  The  notation  /,  g  is  intended  to  be  suggestive,  but  these  filters  are  not,  even  in  a 
discrete  sense,  dual  vectors.  We  also  remark  that  the  specific  implementations  considered 
are  for  a  real  signal  s.  That  is,  we  seek  to  recover  the  signal  from  the  approximation 
s  i?cs^,  that  is,  the  real  part  of  the  DWT'^^  (cf.,  equation  (4.1)).  A  wide  variety  of 
approximate  inverses,  several  of  which  are  variations  on  the  standard  frame  approximation 
(1.16),  are  obtained  simply  by  r^arying  these  filters.  We  shall  find  that  the  DWT~^  can 
achieve  better  results  than  (1-16)  with  comparable  or  even  less  computation.  (A  price  is 
paid,  however,  in  the  forward  transformation  since  we  require  the  undecimated  transform.) 
Moreover,  unlike  (1.16),  all  such  inverses  are  translation  invariant.  That  is,  the  composition 
DWT~^  o  DWT  commutes  with  translations. 


A.  Filter  Constraints 

Clearly,  an  attempt  to  treat  banks  of  arbitrary  filters  would  be  much  too  general  and 
unlikely  to  produce  a  very  useful  theory.  In  this  subsection,  we  propose  some  minimal 
constraints  on  the  forward  and  backward  filters.  At  the  very  least,  the  DWT  should  be 
interpretable  as  a  time-scale  decomposition,  should  be  numerically  .stable,  and  should  have 
a  stable  inverse.  For  example,  it  is  standard  to  require  that  /,  <?  be  a  lowpass-highpass 
filter  pair,  that  /  be  real,  and  that  f„  =  s/2  before  considering  the  forward  filter  bank 
a  DWT  [2].  An  analysis  of  numerical  behavior  is  beyond  the  scope  of  this  work;  however, 
we  do  present  conditions  for  finite  energy  and  boundedness  of  the  DWT  and  DWT~^  as 
the  number  of  stages  M  becomes  infinite.  These  properties  are  important  inasmuch  as  they 
reflect  directly  on  the  numerical  stability  of  the  algorithms  themselves  [4].  Even  though 
one  never  actually  computes  an  infinite  number  of  octaves,  the  poor  numerical  behavior 
of  unbounded  transformations  tends  to  persist  for  finite  M. 

Although  an  appropriate  expression  for  the  energy  of  the  discrete  nonorthogonal  trans¬ 
form  is  somewhat  arbitrary,  we  take  one  which  is  consistent  with  energy  preservation  in 


i.u_  t.: _ 4.1 _ I _  _ j  _ : _ _  _ _ _ 
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require  the  existence  of  Ce  >  0  satisfying  (cf.  equation  (4.26)) 
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r  Zw  10 


We  then  define  the  energy  by 
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This  is  to  be  interpreted  to  mean  that  (lyc^w'”,  s*^)  is  an  energy  density  for  the  DWT 
and  the  total  energy  is  given  by  applying  the  metric  given  by  the  {ML  -t-  l)-diinensional 
diagonal  matrix  with  elements  Although  it  may  seem  strange  that  the 

defiiution  of  E  in  equation  (5.2)  depends  on  the  filters,  the  offending  constant  Cg  is  needed 
to  balance  the  energy  between  /  and  the  Alternatively,  we  could  have  absorbed 

into  the  filter  g  or  used  a  different  mother  wavelet  V/j  ^  ^  Suppose  that  /  is  a 

lowpass  filter  in  the  sense  that  |/^  (aj)|  <  \/2  for  w  ^  0  (recall  that  /^(O)  =  V2)-  Then, 
the  last  term,  which  represents  the  energy  in  the  spectral  interval  [0,  1/2'^j,  goes  to  zero 
as  M  goes  to  infinity  (cf.  Appendix  B).  Thus,  the  constant  represents  a  constant  factor 
in  the  non-DC  energy  . 

Perhaps  the  most  significant  aspect  of  the  definition  (5.2)  is  the  fact  that,  for  M  =  oo, 
the  discrete  transform  is  bounded  and  has  a  bounded  inverse  if  and  only  if  there  exist 
A,  5  >  0  such  that 

Al|s|l^  <  E(w)  <  i5||s||^  .  (5.3) 

We  point  out  that  (5.3)  is  actually  a  discrete  admissibility  condition.  The  following  is  a 
sufficient  condition  for  it  to  be  satisfied  [2] 


0  <  A  < 


<  B  <  oo  for  all  u>  . 


(5.4) 


At  w  =  0,  (5.4)  is  also  a  necessary  condition  (see  Appendix  B).  Hence,  we  must  have 

1 


-JfMV  <  1 


(5.51 


Note  that,  since  /^(O)  =  Yin  ~  condition  (5.4)  requires  Y2n  s'n  — 

that  all  the  o''  be  highpass  filters.  The  DWT  will  be  energy  preserving  (that  is,  E  ~  llslp) 

if  and  only  if  (5,4)  holds  with  A  =  .B  =  1  [2]. 

Finiteness  of  the  DWT~^  as  M  — >  oo  also  imposes  some  constraints.  The  forward 
and  inverse  transforms  include  terms  containing  the  filters  f  and  /  iterated  M  times 
respectively.  We  cannot  permit  the.se  terms  to  grow  without  bounds.  More  precisely,  the 
output  of  the  DWT~^  is  given  by  (cf.,  figure  2b) 


t-l 


i=o  j=o 


M-] 


*=  EH  [(DV)  *  ]  (D*y)  *  v/‘  +  n  [(»'/) 


>=0 


(5.6) 


For  i  =  0,  the  first  term  on  the  right  is  understood  to  mean  y  >*- w''.  Also,  equivalent  to  the 
condition  that  rj}{t)  be  real  or  analytic,  we  require  that  either  both  g  and  g  are  real  or  one 

^This  would  make  some  .sense,  since  it  makes  the  role  of  /  more  explicit  in  the  relationship  of  the 
continuous  and  discrete  transforms.  Each  filterbank  has  a  “true”  mother  wavelet  associated  to  it, which 
only  approximates  the  original  See  [2],  [lOj. 


of  them  is  analytic  (that  is  has  zero  z-transform  for  w  <  0).  If  the  inverse  is  to  be  bounded 
in  the  metric  of  equation  (5.2),  then  uniformly  bounded  inputs  w',  must  produce 
uniformly  bounded  outputs.  Taking  z-transforras,  we  find  this  entails  that  J  |s"  (a;)|V2- 

uniformly  bounded  in  M  implies  that  /  rijLo*  l/*(2'’'^))P  |s^^(a’)p  is  uniformly  bounded. 
This,  in  turn,  implies 

l/z(0)i  <  ^  .  (5.7) 

Similarly,  if  the  composition  of  the  forward  and  inverse  transforms  is  to  be  bounded  for 
arbitrju'y  M,  then  we  must  have 

l/,(0)/,(0)i  <  1  .  (5.8) 

Note  that,  since  we  require  that  /^(O)  =  \/2,  (5.8)  implies  (5.7).  Details  and  conditions 
on  the  inverse  filters  similar  to  (5.4)  may  be  found  in  Appendix  B. 

Of  course,  it  would  be  desirable  that  the  DWT~^  be  an  exact  left  inverse.  As  described 
in  the  introduction,  this  will  be  true  if  and  only  if  equation  (1.14)  holds.  For  multiple  voices 
and  complex  wavelets,  this  takes  the  form 

Re'^g'’  *  g"  +  J*f  =  6  .  (5.9) 


If  /  and  g  are  filters  corresponding  to  orthonormal  wavelets,  then  (5.9)  will  bold  for 
a  single  voice  with  /  =  /^  and  g  =  g^ .  Unfortunately,  in  the  case  of  Morlet  wavelets, 
backward  filter.s  which  satisfy  (5.9)  are  generally  much  too  long  for  most  applications, 
In  contrast,  even  in  the  nonorthogonal  case,  we  shall  find  that  the  adjoint  filters  and  g^ 
provide  useful  approximations.  Other  backward  filters  of  interest  are  /  =  c6  and  g  =  c6. 
These  filters,  which  are  motivated  by  the  continuous  inverses  (4.1)  and  (4.15)  will  prove 
extremely  effective.  The  constraint  (5.7)  takes  the  form 

/  =  c implies  c  <  1/2  ,  (5.10) 


and 


f  =  c  S  imnli.es  e  <  1  /  v/2 
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B.  Adjoint  DWT 

We  pause,  briefly,  to  show  that  the  filter.s  f  ■—  and  g  =  g^  yield  the  adjoint  DWT. 
Let  H*  and  be  the  linear  transformations 


H‘  ^  J]  [(D7)  * 

J-O 


{D'g)  * 


M-i  _ 

and  L"  =  n  [(»’?)* 

j=u 


(5.12) 
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From  figure  2a,  the  forward  transform  DWTfg  is  given  by  w'  =  H's  arid  s"  =  L*^s.  For 
the  Euclidean  metric  we  have 

<w,  DWTfg  s>  =  wj,  |H^]“  +  w;;  [L«s>^]„ 


=:  53 

i.n 

=  5]  <H' V,  s>  +  s> 


=  <DWT-, \,w,s>  , 


(5.13) 


where  the  last  line  follows  from  (5.6)  with  f  —  and  g  =  .  This  implies  that  DWT.^^  ^ 

is  identical  to  (DWTfg)' .  Our  algorithm  actually  uses  /  /2  which  also  satisfies  (5.8), 

It  is  easy  to  see  that  this  yields  the  ad.joint  under  the  energy  metric  (5.2). 

C.  Implementations 

In  addition  to  the  baffling  number  of  inverses  a\'ailablc  for  consideration,  one  of  the 
most  confusing  aspects  of  this  subsection  is  a  repeated  shifting  between  decimated  and 
undecimated  wavelets.  The  gist  of  these  machinations  is  that,  although  our  wavelet  trans¬ 
forms  as  well  as  the  inverse  filterbanks  will  always  be  undecimated,  the  interpretation  of 
and  distinctions  between  the  various  inverses  are  more  properly  made  by  an  appeal  to  the 
decimated  case.  Much  of  the  misery  disappears  if  one  keeps  in  mind  that  [w*']„  ;=  [w"’]2in; 
that  iS,  the  two  tj'pes  of  Wavelet  tiaiisfoiiii  are  simply  uiueient  sainpling.s  of  the  same 
function  and  that  they  coincide  at  their  common  points.  We  also  remind  the  reader  that 
v/hen  y  is  complex,  s  corresponds  to  the  real  part  of  s,  that  is,  Res.  For  notational 
simplicity,  in  this  subsection,  we  present  our  derivations  using  only  a  single  voice  and 
replacing  the  constant  2ln2/ C^j,  by  c.  The  implementations  of  three  inverses  will  be  dis¬ 
cussed:  the  double-integral  type  resembling  (4.3a),  the  adjoint  DWT,  and  a  single-integral 
type  corresponding  to  (4.17).  The  first  two  bear  an  intimate  relationship  to  the  frame 
approximation. 

C.l.  Double-integralJijpe 

We  bevin  with  f  =  /i  I ^/2  ;^.nr1  n  =  AHVirnviftfincr  flip  t.prm  nf  Kv  Hn  wp 

-  -  - /  V  —  - - - - cp - -  - - -  V'"  ’'/  's/  “■^7 

have 


=  <=  E  4:  [(D's')  •  w]„  +  DC 

i=0  ^ 

M-1 

=  ^  E  E^E^X.n  +  DC 

i=0  n 


i-O  ^  ^  n 
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(5.14) 


M-l  .. 

i=D  n  ^  ^ 

M-l 

=  C  ^  5]]wi,„V»i,n(0)  +  DC  . 

i=0  n 


Except  for  the  DC  term,  this  coincides  with  the  standard  approximation  (1.16)  at  time 
t  =  0. 

We  wish  to  extend  this  formula  to  arbitrary  times  t  m.  Let  w*(m)  =  w*(T_ms) 

and  Wi(m)  =  Wj(T_,„s)  be  the  wavelet  transforms  obtained  by  first  translating  the 
signal  by  m  points  (cf.  (2.9)).  The  idea  is  that  these  transforms  are  computed  at  (i.e., 
centered  about)  time  m.  In  this  context,  the  standard  approximation 


Af-l  M-l 

Sm  =  s{m)  JJi  Rec  ^  ^ Wi,„V>,,„(m)  =  ijec  ^  ^ Wj,„(0)i/>.,„(m)  (5.15) 

1=0  n  1=0  n 

utilizes  wavelet  coefficients  “computed”  at  time  0.  For  m  —  0,  equation  (5.15)  gives 
So  ~  jFie  X- ^un(^)V't,n(0).  Also,  it  is  apparent  from  the  form  of  (5.15)  that  for  m  ^  0, 
the  decimated  wavelets  extrapolate  s(t)  from  time  t  =  0  to  time  t  =  Since 

■Wi,n+2^(0)  —  ''''^«,n(2'^“‘),  the  extrapolation  is  never  longer  than  2'^  —  1.  Alternatively,  we 
could  translate  the  signal  by  m,  compute  Wi,„  on  the  translated  signal  T-n^s  (which  by 
definition  yields  wgn(’T^'))>  evaluate  the  expansion  at  time  0.  This  yields 


M-l 


~  [T-ms]o  ^  Re  c  EE  w,,„(m)t/>i,„(0) 


(5.16) 


t=0  n 


thereby  avoiding  any  extrapolation.  Note  that  (5.15)  and  (5.16)  represent  two  quite  dif¬ 
ferent  approximations.  This  is  because  the  decimated  DWT  is  not  translation  invariant; 
i.e.,  Wj^„(m)  ^  Wi_„+,n(0).  In  fact,  we  would  expect  (5.16)  to  be  a  better  approximation 
since  it  does  not  extrapolate.  (The  tradeoff  is  tha.t  we  must  compute  Wi(')  at  every  point, 
which  is  equivalent  to  computing  the  undecimated  DWT.) 

On  otliGr  hsiiid,  tlic  und^cimstcd  trcinsforiri  stncl  its  inverse  2.rG  tim^^  invcirisint. 
This  implies  that  =  [T-^iSjo.  Substituting  this  into  (5.14)  yields 


M-l 

Srr,  =  [T_,„s]o  =  C  ^  ^  Wi,„(m)  V'<,n(0)  +  DC  ,  (5.17) 

«  =  0  n 


^®V/e  use  the  word  extrapolate  to  emphasize  that  the  contribution  of  a  single  wqn  i-^  exti  ajiolated  via 
V’i.n.  On  the  other  hand,  since  in  a  sense  we  are  spanning  the  gap  of  2'  points  between  Wj-^n  and  w,  , 
which  in  turn  depend  on  the  entire  signal  s,  wc  are  interpolating.  This,  of  course,  is  a  much  niore  stable 
operation  than  strict  extrapolatien. 


C«U 
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» 
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which,  except  for  the  DC  term,  coincides  with  (5.16).  In.  summary,  we  have  considered  two 
approximations  to  s[m):  (i)  the  decimated  transform  computed  at  time  0  and  employed  in 
the  corresponding  wavelet  expansion  to  approximate  the  signal  at  time  m  (equation  (5.15)), 
and  (ii)  tie  decimated  transform  computed  at  time  m  used  in  the  wavelet  expansion  to 
approximate  s{m)  —  [T-ms]o  at  the  same  time  as  the  computation  (equation  (5.16)).  As 
the  numbe.i  of  scales  M  goes  to  infinity,  the  latter  formula  becomes  identical  to  (5.17),  that 
is  to  the  DWT~^  (i-e-i  undecimated  inverse  filter  bank)  with  /  =  /v^i  and  g  =  gK  In 
fact,  (5.17)  represents  an  improvement  over  (5.16)  since  the  DC  term  replaces  the  energy 
lost  in  truncating  the  expansion  at  i  =  M  —  1. 

C.2.  Adjoint 

We  already  know,  from  (5.13),  that  the  DWT~^  with  f  =  /V^  and  5  =  is  the 
adjoint  to  the  DWT.  It  also  has  a  simple  interpretation  in  terms  of  the  frame  approximation 
but  which  is  a  bit  difficult  to  demonstrate.  It  is  an  average  of  approximatio.ns  of  the  form 
(5.15).  More  precisely,  let  *  D  be  the  operator  which  acts  on  vectors  x  by  (/^  *  D)x  = 
♦  (Dx).  As  discussed  in  [2],  when  /  is  an  a  trous  filter  this  operator  is  a  dyadic 
interpolator.  That  is,  [/^  +  (D^^)]m  «  so  that 

[(/t*D)Vi,„  •  (5.18) 


Furthermore,  it  can  be  shown  (Appendix  E)  that  the  DWT  ’  with  /  ~  /2  and  g  =  cg^ 

yields 

M-l  2'-l 


u  ft 

1=0  r=0 


n\«  ^  I 


Substituting  (5.18)  in  (5.19),  we  have 


M-l  ^  2*-l 

Wi,„(r)  -  r)  +  DC  (5.20) 

i=0  n  r=0 

where  we  have  used  (E.7).  One  rna.y  rewrite  this  as  c  '  ^^n+r(’^0  reflect- 

imr  the  annroximation  i/)l  fe  ci/;!  of  f4.4al.  We  ma.v  consider  the  set  of  —  r)  — 

W  AA  f  \  -  /  ^  -  -  --  -  -  __  > 

^  ^  ~’^k  with  k  =  2'n  +  r  of  (4.11a).  On  the  other  hand,  the  dec¬ 

imated  interpretation  is  somewhat  more  enlightening.  The  decimated  transform  at  oc¬ 
tave  i  is  invariant  under  translations  of  2';  i.e.,  Wi,„(2*)  =  w,_nqi(0).  For  each  value  of 
r  =  0, . . . ,  2*  —  1,  we  can  use  (5.15)  to  estimate  by  extrapolating  from  time  r  to  time 
m,  that  is,  Sm  ~  7?e  Equation  (5.20)  is  simply  an  average 

of  these  extrapolations.  Since  the  longer  the  extrapolation,  the  less  likely  we  are  to  have 
a  good  approximation,  the  quality  of  this  result  should  lie  somewhere  in  between  that  of 
(5.15)  and  (5.17).  Hov/ever,  since  the  double-integral  type  (5.17)  not  only  provides  a  better 
approximation,  but  also  requires  much  less  computation  (a  scalar  multiplication  for  6l\/2. 


versus  a  convolution  for  the  filter  /V2),  there  seems  to  be  little  reason  to  use  the  adjoint 
inverse. 


C.3.  Single-integral  type 

Finally,  we  examine  g  =  c6  with  /  =  6/\/2.  This  results  in  an  inverse  analogous  to 
the  single  integration  approximation  (4.17).  More  precisely,  substituting  ~  ^/\/2  and 
=  cS  into  (5.6),  we  have 


s 


m 


M-1 


=  E 

»=o 


(5.21) 


Essentially  the  same  formula  holds  for  decimated  wavelets;  i.e.,  with  replaced  by  Wi^m- 
However,  in  that  ceise  not  all  octaves  are  available  at  each  time.  More  precisely,  if  we  set 
m  =  2®r  +  2”n  where  2‘^r  <  2’'^  and  Q  <  M  is  the  largest  integer  satisfying  such  an 
expression,  we  have 


i=0 


(5.22) 


If  one  omits  the  DC  term,  this  formula  will  provide  very  poor  results,  since  at  those  points 
for  which  Q  is  small,  up  to  one-half  the  energy  may  be  omitted.  Since  it  is  common 
practice  to  use  only  the  wavelet  coefficients  Wj,  this  may  explain  the  reticence  of  many 
people  to  use  the  single  integration  approximation.  In  contrast,  (5.21)  yields  excellent 
results,  often  better  than  (6.17)  or  (5.2u)  (cf.  Section  VI).  Moreover,  this  DWT~^  does 
not  increase  the  temporal  support  beyond  that  of  the  wavelet  coefficients,  a  property  which 
is  extremely  helpful  for  filtering  in  the  wavelet  domain.  These  features,  coupled  with  a 
high  computational  efficiency,  make  this  inverse  our  choice  for  most  applications. 


D.  Normalizing  g^ 

In  Section  IV,  we  related  particular  realizations  of  the  DWT~^  to  discrete  approx¬ 
imations  to  the  continuous  inverse.  This  correspondence  was  determined  up  to  the  DC 
term  and  a  constant,  i.e.,  s  ps  cs.  In  theory,  this  gives  us  two  constants  to  adjust;  one 


u  f  «•§•*>*» 

Afc*  ^  VX.'A  AAA 


AVA  VAA^  X^VAA^A  AAA  VAA  UA^  A  A  AA^  VAAl.^ 


nni,  Aa**  >3 1  1  CiCi  t  1 

A.  A/.V^  ^VvAAV^A  AA/A  U  A  WA  A 


how  to  pick  these  constants  (or  for  that  matter  more  general  backward  filters)  is  complex, 
application  dependent,  and  beyond  the  scope  of  this  study.  Here,  we  restrict  ourselves  to 
just  a  few  ca.ses.  First,  to  avoid  a  dependence  on  M,  we  determine  the  constants  from 
the  consideration  of  a  generic  single  stage  of  DWT~^  o  DWT.  NJe  shall  remark  on  two 
possibilities,  both  involving  approximations  to  (5.9).  The  most  natural  is  to  minimize  the 
error  in  the  least  squares  sense.  That  is. 


min  llciie 

CfC'  ■  * 

V 


(5.23) 
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where  the  form  of  the  backward  filters  is  given  by  one  of  the  implementations  in  the  last 
subsection  and  the  subscript  U  indicates  the  renormalizations  y  =  cg^  and  f  —  c'fjj. 
We  observe  that  if  the  zero‘d  component  of  (5.9)  holds  exactly  ,  then  the  substitution 
9u  =  and  fy  =  j2  yields  c  J2v  +  (c' /2)  [f^  *  /]o  =  1.  That 

isi  ^  Jlv  T”*'  +  (‘^72)11/11'*  =  1-  This  is  simply  (4.26)  with  c  =  (2\n2,)/(C^L)  and 

c'  =  1. 

A  second  approach  is  to  require  that  the  energy  of  the  reconstruction  match  that 
of  the  signal.  To  do  this  in  a  signal  independent  manner  is  not  generally  possible.  As 
a  compromise  we  apply  this  re  yairement  to  a  specific  signal,  namely  an  impulse,  which 
has  a  fiat  spectrum.  Setting  the  energy  of  the  impulse  response  of  a  single  stage  of  the 
DWT~^  o  DWT  equal  to  one  yields  j|ciEe  f  v  *  fW'^  ~  Note  that, 

because  of  the  cross- terms,  this  differs  from  preserving  the  energy  of  the  wavelet  transform 
(cf.  equation  (5.1)).  Rather,  it  attempts  to  take  into  account  the  nonorthogonality  of 
f  *  f  and  g  *  g-  Contrary  to  intuition,  informal  trials  indicate  that  it  is  a  better  choice 
than  (5.23);  i.e.,  ||i2es  —  s|p  is  smaller.  It  cannot  be  better,  of  course,  for  the  case  in 
which  M  =  0  and  the  signal  is  impulsive,  but  it  does  seem  comparable  and  provides  an 
improvement  for  M  >  0.  Finally,  we  have  found  that  taking  c'  =  1  seems  to  work  as  well 
or  better  than  any  other  A^alue,  so  we  shall  accept  that  simplification.  With  the  constant 
c  included  in  the  g"'s,  this  condition  then  becomes 


+  7*/f  -  1  •  (5.24) 


The  author  postulates  that  other  choices  will  at  best  provide  a  marginal  improvement 
and  that  if  more  accuracy  is  desired,  it  is  much  more  efficient  to  iterate  the  suggested 
implementations  by  using  Algorithm  3.1  than  to  look  for  other  normalizations  of  /  and  g. 
Another  alternative,  which  works  very  well  in  practice  and  which  will  also  be  evaluated 
in  Section  VI,  is  to  determine  c  from  the  continuous  approximations  (4.3)  and  (4.17)  by 
which  5  =  {2\n2)/C,p)g^  and  g  =  (2  In 2)/C'i^)^  respectively. 

Let  us  combine  these  remarks  with  the  results  of  subsection  V.C  to  summarize  the 
various  DWT~^'s  under  consideration.  We  also  shall  take  advantage  of  this  opportunity 
to  add  voices  to  our  formulation  wherever  they  were  previously  omitted.  There  are  two 
choices  of  g'’  under  consideration.  The  double  integral  type  is  given  by 


^  ('/.’’It.  7  j:/,/9 


yV 


\9  /  y  ^ 


with  C'2  determined  by 


L-i 


u=0  ^ 


(5.256) 


Alternatively,  by  analogy  to  the  continuous  case,  we  may  set 

21n2 

C2  W 


L  Cff, 


(5.25c) 


The  voice  normalization  7  *'  in  (5.25a)  comes  from  generalizing  (5.14)  along  the  lines  of 
(4.10b). 

The  single  integral  type  (cf.  (4.18))  is  given  by 


with  Cl  determined  by 


L-\ 


=  1 


or  by  analogy  to  the  continuous  case 


21n2 


Cl 


LC 


(5.26a) 


(5.266) 


(5.26c) 


One  might,  quite  justifiably,  wonder  about  the  absence  of  1/2*  in  (5. 25)  or  (5.26) 
when  comparing  them  with  (4.11)  or  (4.18).  In  the  case  of  (5.25),  it  is  (4.11b)  not  (4.11a) 
which  is  relevant.  More  precisely,  the  reader  is  reminded  that,  although  the  input  to  the 
DWT''^  is  the  undecimated  tra.nsform  wjj",  its  output  (5.17)  corresponds  to  a  time-shifted 
decimated  inverse;  that  is,  to  (4.11b)  which  has  no  1/2*.  The  filter  /  =  6/ ^'2  does  supply 
a  factor  of  However,  it  is  absorbed  into  the  normalization  of  „  as  detailed  in 

equation  (5.14).  This  contrasts  with  the  adjoint  DWT~^  where  f  —  12  results  in  an 

additional  factor  of  1/2  (/o/2  =  2'^^^  versus  2~’/^)  which  actually  appears  in  (5.20)  so 
that  the  adjoint  DWT~^  ’s  output  corresponds  to  (4.11a).  Finally,  in  the  DWT~^  specified 
by  (5.26),  /  =  0/V2  is  repeated  i  times  producing  a  factor  of  1  j^'2'  which  appears  in  (5.21) 
and  corresponds  to  (4.17).  A  simple  way  to  normalize  /  is  to  take  equality  in  (5.7);  that 

is,  En7n  =  l/N/2- 


E.  Summary 

The  following  table  summarizes  the  DWT~^  implementations  treated  in  this  section. 


Table  5.1  DWT-^ 


DWT-^ 

7 

a 

s 

Interpretation 

Adjoint  (Ad) 

Double  Integral  (DI) 

Single  Integral  (SI) 

2'' 

sff 

J 

(5.20) 

(5.17) 

(5.21) 

adjoint  of  (1.11) 

approx  (4.3b)  to  (4.1) 

approx  (4.17)  to  (4.15) 

i 


i 
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Remarks 

(a)  The  DWT~^  is  computed  using  undecimated  wavelets  w*  according  to  figure  2. 

(b)  The  relationship  between  in  the  referenced  equations  and  the  filters  g'’  is  given 
by  gn  =  ^(-n)  = 

(c)  In  cases  Ad  and  DI,  it  is  assumed  that  g  is  real  or  that  it  is  aixalytic;  i.e.,  ^  0 

for  w  <  0. 

(d)  In  case  SI,  it  is  assumed  that  g  is  analytic  or  that  =  g  (i.e.,  is  real).  Note 

that  neither  condition  holds  for  the  Daubechies  orthonormal  wavelets. 

(e)  The  definitions  of  the  constants  C2  and  cj  are  found  in  (5.25)  and  (5.26)  respectively. 

(f)  The  factors  1/2  and  C2/7*’  appearing  in  the  filters  of  the  adjoint  imply  that  it  is  the 
adjoint  with  respect  to  the  energy  metric.  See  equation  (5.2). 

(g)  In  the  case  of  an  a  trous  implementation  of  Morlet  wavelets,  /  =  /^  is  real  symmetric 
and  {g'^'y  =  g  is  complex  Hermitian.  In  that  case,  the  f’s  may  be  omitted  in  table  5.1. 

(h)  The  equations  cited  under  s  involve  decimated  wavelet  coefficients  even  though  s  is 
equal  to  the  output  of  the  undecimated  inverse.  See  the  introduction  to  subsection 
V.C. 

(i)  The  abbreviations,  which  will  also  be  used  in  Section  VI,  are  Ad=Adjoint,  DI=DoubIe 
Integration,  SI=Single  Integration. 


VI.  COMPARATIVE  PERFORMANCE  OF  INVERSES 


In  this  section,  we  compare  the  three  DWT~^  filter  banks  featured  in  table  5.1.  To 
this  we  add  the  standard  frame  approximation  (5.15).  Tiie  task  is  considerably  more 
complex  than  might  at  first  appear.  The  performance  of  these  algorithms  is  bound  to  be 
both  wavelet  (including  bandwidth,  number  of  octaves  and  number  of  voices)  and  signal 
dependent.  In  addition,  there  are  variations  on  the  algorithms  which  we  would  like  to 
investigate,  such  as  the  choice  of  cj  and  C2  or  the  inclusion  of  the  DC  term.  Also,  one 
may  apply  Algorithm  3.1  to  these  DTVT'^’s  resulting  in  true  left  inverses.  The  quality 
of  the  particular  DWT~^  then  becomes  a  question  of  the  speed  of  convergence  of  the 
iterations.  One  perfectly  legitimate  question,  that  of  the  relative  behavior  of  these  inverses 
after  filtering  in  the  wavelet  domain,  is  beyond  the  scope  of  this  paper. 

In  order  to  put  some  reasonable  bounds  on  what  will  be  an  essentially  empirical  study, 
we  restrict  ourselves  to  Morlet  wavelets  of  the  form 


and  consider  only  two  signals,  a  discrete  impulse  and  a  sinusoid.  They  are  the  prototype 
examples  of  two  major  classes  of  signals,  narrowband  and  broadband.  The  former  outputs 
the  transfer  function  of  the  linear  system  DWT~^  o  DWT  which  may  then  be  com2:)ared  to 
that  of  a  true  left  inverse,  in  this  case  a  delta  function.  Except  to  examine  the  single  stage 
behavior,  we  shall  fix  the  number  of  octaves  at  three,  that  is,  i  =  0,  1,  2.  The  frequency 
bands  of  these  octaves  are  [0.25,  0.5],  [0.25,  0.125],  and  [0.125, 0.0625].  The  sinusoid  which 
we  shall  use  is  sin(O.lt),  sampled  at  At  =  1.  Consequently,  its  energy  is  concentrated  in 


the  last  octave.  Other  sinusoids  yield  siin'l...!  results;  however,  it  should  be  mentioned  that 


the  placement  of  the  sinusoid  frequency  relative  to  the  centers  of  the  voices  y"  does  make 
a  difference.  Also,  to  avoid  the  effects  of  the  onset  and  termination  of  the  input  (which 
resemble  impulses  and  thus  are  broadband),  the  sinusoidal  signal  is  generated  from  0  to  T, 
but  only  data  in  the  interval  [T/4,  3T/4]  are  used  for  comparing  the  signal  to  the  wavelet 


inverse. 


An  extremely  important  aspect  of  this  study  is  the  choice  of  error  criterion.  The  most 
obvious  is  to  take  a  function  of  the  difference  s  —  m  Although  this  works  quite  well  in  the 
case  of  the  imi^ulse,  problems  occur  for  the  sinusoid.  The  algorithms  produce  an  accurate 
sinusoidal  inverse,  but  they  do  not  preserve  the  amplitude;  that  is,  they  reproduce  the 
signal  up  to  a  multiplicative  factor.  For  example,  for  a  sinusoidal  signal  of  unit  norm  and 
one  typical  DWT~^  o  DWT,  the  root  mean  square  error  for  the  raw  inverse  was  0-125, 
but  after  the  inverse  v/as  normalized  to  have  unit  norm,  the  error  was  only  0.007.  Under 
the  same  circumstances,  the  rms  errors  for  the  impulse  were  0.33  and  0.34  respectively. 
Consequently,  in  order  get  a  measure  of  the  quality  of  the  inverse  exclusive  of  the  amplitude 
factor,  we  shall  base  our  criteria  on  normalized  signals. 
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Our  error  function  shall  be  the  root  mean  square  error 


^rms  — 


S  —  S' 


(6.3) 


where  s'  =  (||sl|/||s||)?.  Another  possible  choice,  for  example,  would  be  the  maximum 
absolute  error  Cmax  —  max^  An-  For  convenience  and  without  loss  of  generality,  the  signal 
will  always  be  chosen  such  that  ||s|j  —  1.  By  also  displaying  the  output  energy  l|s|p,  which 
for  a  true  left  inverse  should  be  1.0,  we  keep  an  eye  on  the  amplitude  mismatch.  We 
emphasize,  however,  that  the  discrepancy  in  the  norms  ||s||  and  ||s[{  is  signal  dependent 
and  simply  renormalizing  the  inverse  does  not  solve  the  problem.  Clearly,  when  the  signal 
consists  of  a  mixture  of  narrowband  and  broadband  components,  the  amplitude  distortion 
introduced  by  the  DWT~^  is  more  than  a  simple  multiplicative  factor. 

In  order  to  dispel  the  visual  vacuum,  we  begin,  in  figure  3,  with  plots  of  three  typical 
inverses.  All  three  reproduce  the  impulse,  but  with  noticeable  side  lobes.  The  large  dip 
at  zero  for  the  frame  approximation  (Pr)  is  due  to  the  lack  of  a  DC  term.  The  only  other 
obvious  qualitative  difference  is  the  smaller  support  of  the  single  integration  (SI)  formula. 
A  great  deal  more  may  be  said  from  an  examination  of  table  6.1.  We  are  immediately 
struck  by  the  high  quality  of  the  inverses  for  narrowband  signals.  Perhaps  this  accounts 
for  the  success  of  the  standard  approximation  in  many  applications.  On  the  other  hand, 
the  behavior  for  an  impulsive  signal  are  often  inadequate.  This  implies  that  iteration  may 
well  be  necessary  and  consequently  should  be  an  important  consideration  in  evaluating  the 
various  inverses.  The  four  inverses  appearing  in  table  6.1  are  presented  in  approximate 
order  of  quality.  The  SI  inverse  has  the  smallest  error  and  is  by  far  the  most  efficient  to 
implement.  Although  only  a  few  parameter  values  are  presented  here,  these  inverses  seem 
to  exhibit  a  similar  relative  behavior  over  a  large  range  of  values  of  /?  and  of  L. 


Figure  3.  The  transfer  functions  {DWT  o  DWT"^)  of  three  selected  inverses.  The  signal  i.";  an  impulse, 
and  the  forward  transform  is  a  Morlet  wavelet  with  P  =  0.125,  four  octaves,  and  8  voices.  The  inverses 
shown  are  (a)  single  integration,  SI;  (b)  double  integration,  DI,  (c)  and  the  frame  approximation,  Fr. 


TABLE  6.1  Comparison  of  Four  Inverses  (3  octaves). 


impulse 

sine 

DWT-^ 

p 

L 

RMS 

ENERGY 

RMS 

ENERGY 

SI 

0.5 

4 

0.31 

1.00 

0.0024 

1.19 

0.125 

8 

0.22 

1.06 

0.039 

0,55 

DI 

0.5 

4 

0.34 

0.88 

0.0069 

0.77 

0.125 

8 

0.26 

0.99 

0.044 

0,49 

FT 

0.5 

4 

0.50 

0.65 

0.0072 

0.67 

0.125 

8 

0.39 

0.87 

0.050 

0,60 

Ad 

0.5 

4 

0.32 

0.96 

0.00S3 

0.77 

0.125 

8 

0.54 

0.78 

0.067 

0.31 

Figure  4  illustrates  the  Neumann  iterations  (cf.  Algorithm  3.1)  for  a  DWT  of  type  SI, 
and  table  6.2  compares  the  results  of  such  iterations  for  the  four  major  types  of  inverses 
when  the  original  signal  is  an  impulse.  Inasmuch  as  the  energy  rarely  increased  under 
DWT  o  DWT~^  and  always  less  than  a  factor  of  2,  /i  was  nonrinally  set  to  0.5.  In  the 
absence  of  other  criteria,  SI  seems  to  be  by  far  the  best  choice.  The  incredibly  slow 
convergence  for  Fr  is  due  to  the  lack  of  a  DC  term.  (Optimizing  /u  made  little  difference. 
The  iterations  of  Fr  do  not  converge  much  faster  for  larger  p  and  diverge  in  the  region 
of  /i  =  3.)  In  fact,  the  same  difficulties  occur  with  the  other  inverses  if  the  DC  term  is 
omitted.  Without  it,  one  must  include  a  sufficient  number  of  octaves  to  encompass  all  of 
the  signal’s  energy  whether  or  not  a  fine  scale  analysis  of  the  lower  frequency  energy  is  of 
interest.  For  an  impulse,  whose  energy  covers  the  entire  spectrum,  an  inverse  without  the 
DC  term  is  particularly  poor.  Note  that  a  DC  term  may  be  added  to  the  Fr  inverse  by 
constructing  the  scale  function  and  including  the  term 


1  ,  /  m  \ 


(6.4) 


in  equation  (5.15).  In  view  of  the  other  inverses  available,  this  does  not  seem  to  be  worth 
the  effort.  It  should  be  mentioned,  however,  that  the  frame  expansion  does  provide  a 
continuous  interpolation  of  the  signal  (i.e.,  replace  m  by  t),  a  property  not  shared  by 
the  DWT~^ .  Finally,  we  remark  that  the  difference  between  the  expressions  (5.25b)  and 
(5.25c)  for  C2  or  between  (5.26b)  and  (5.26c)  for  ci  seems  to  be  insignificant. 
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(a)  (b) 

Figure  4.  The  result  of  performing  Neumann  iterations  on  the  inverse  SI  after  (a)  5  iterations  and  (b) 
10  iterations.  The  forward  transform  had  /S  =  0.125  and  T  =  8;  the  iteration  factor  was  p  =  0.5. 

TABLE  6.2  RMS  errors  of  iterated  DWT~^  with  Morlet  wavelets  and 
implulsive  signal;  fi  =  0.5,  3  octaves,  L  =  4,  fi  =  0.5. 


NUMBER  OP  ITERATIONS 

DWT-'^ 

1 

5 

10 

25 

50 

100 

SI 

0.264 

0.120 

0.048 

0.0038 

0.000067 

2.9  10-^ 

DI 

0.308 

0.190 

0.115 

0.033 

0.0050 

0.00016 

Pi- 

0.461 

0.371 

0.322 

0.284 

0.273 

0.264 

Ad 

0.283 

0.162 

0.091 

0.020 

0.0021 

0,00012 

The  number  of  voices,  which  is  really  a  parameter  of  the  forward  transform,  is  par¬ 
ticularly  relevant  to  the  quality  of  the  inverse.  The  optimal  number  of  voices  (for  a  given 
0)  tends  to  be  that  given  by  equation  (4.29)  with  equality,  that  is,  approximately  1//?. 
Too  few  voices  greatly  degrades  the  results.  Without  iteration,  additional  voices  make 
little  difference.  However,  a  degradation  in  the  speed  of  convergence  of  the  Neumann  it¬ 
erations  does  seem  to  take  place  with  additional  voices.  For  example,  let  /?  =  0.25,  for 
which  the  suggested  number  of  voices  is  L  =  3.  Then,  the  mse  errors  of  SI  at  iteration 
10  for  L  =  1,  2,  3,  8  are,  respectively,  0.25,  0.018,  0.0032,  and  0.028.  This  phenomenon 
is  somewhat  surprising  inasmuch  as  adding  voices  tends  to  make  a  tighter  frame  [3],  [4]. 
However,  ultimately,  convergence  hinges  on  the  specific  eigenvalues  of  the  functional  itera¬ 
tion,  which,  in  turn,  depend  on  the  filter  /  and  on  u  as  well  as  the  number  of  voices.  The 
value  used,  p  =  0.5,  is  well  below  the  point  of  instability,  which  for  the  above  range  of 
parameters  occurs  between  p  =  1.0  and  ,u  =  3.0.  The  best  convergence  seems  to  occur  near 
instability  (although  there  is  a  short  interval  just  before  that  point  where  it  slows  down). 
Further  study  would  be  needed  to  fully  explain  the  situation.  Note  that  the  minimal  values 
satisfying  (4.29)  for  /3  =  0.5  and  /?  =  0.125  are  L  =  1  and  L  =  G  respectively.  However, 
since  for  visual  presentation  it  is  a  good  idea  to  use  at  least  four  voices,  and  since  it  is 
convenient  to  pick  L  »  1//3,  we  have  used  L  ~  4  and  8  in  the  above  presentation. 
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VII.  CONCLUSION 


The  standard  inversion  procedure  for  discrete  nonorthogonal  transforms  is  a  finite 
expansion  in  terms  of  the  analyzing  wavelet.  While  this  so-called  wavelet  series  or  frame 
approximation  works  quite  well  for  many  signals,  it  fails  to  achieve  good  accuracy  or 
requires  an  excessive  number  of  scales  for  others.  To  remedy  the  situation,  we  have  pro¬ 
posed  an  inverse  filter  bank  DWT~^  as  a  prototype  inverse  discrete  wavelet  transform. 
It  provides  a,  unifying  framework  under  which  the  various  (approximate)  left  inverses  are 
obtained  by  varying  the  filters.  For  example,  it  was  shown  that  the  adjoint  transformation 
of  the  undecimated  discrete  wavelet  transform  is  computed  by  using  the  adjoints  of  the 
filters  from  the  forward  transform.  More  generally,  within  this  context,  we  have  exam¬ 
ined  the  properties  and  interrelationships  of  dual  wavelets,  metrics,  pseudo-inverses,  and 
discretizations  of  continuous  inverses.  In  particular,  we  have  found  that  the  absence  of 
decimation  results  in  dual  wavelets  that  are  time  invariant,  and  that  the  standard  frame 
approximation  may  be  interpreted  as  an  interpolation  of  a  DWT~^ .  Finally,  in  addition 
to  the  adjoint,  we  proposed  two  other  discrete  left  inverses,  which  use  filters  based  on 
continuous  integral  formulae. 

Two  major  differences  stand  out  between  the  standard  frame  approximation  and  the 
inverses  discussed  in  this  paper.  They  are  the  use  of  the  undecimated  wavelet  transform 
as  input  to  the  DWT~^  and  the  inclusion  of  a  DC  (low-frequency)  term  which  is  lacking 
in  the  fr2ime  approximation.  Avoiding  decimation  does  provide  redundancy  and,  hence, 
a  potentially  better  approximation.  However,  our  preference  for  not  decimating  actually 
stems  from  the  forward  transform.  A  desire  to  achieve  finer  output  resolution  coupled 
with  time  invariance  makes  the  undecimated  DW^T  the  natural  settinp'  for  nonorthogonal 
wavelets.  Given  this  situation,  it  makes  sense  to  use  an  inverse  that  exploits  all  the  numeric 
information  available.  The  lack  of  a  DC  term  in  the  frame  approximation  is  a  more  serious 
issue.  It  leads  to  poor  accuracy  in  the  case  of  broadband  signals  and  even  implies  that  the 
composition  of  the  forward  and  inverse  transforms  is  singular  so  that  iteration  need  not 
converge  to  a  true  left  inverse. 

The  best  among  the  inverses  studied  seems  to  be  the  single-integral  inverse.  This 
inverse  is  equivalent  to  simply  summing  the  weighted  wavelet  coefficients  across  scale.  It  is 
by  far  the  most  efficient  computationally,  is  generally  the  most  accurate,  converges  much 
more  rapidly  than  others  under  iteration,  and  does  not  increase  the  support  of  the  DWT. 
The  latter  property  is  quite  useful,  inasmuch  as  it  implies  that  filtering  in  the  wavelet 
domain  is  a  quasi-local  operation  in  time.  That  is,  the  value  of  the  DWT~^  at  a  given 
time  only  depends  on  the  wavelet  coefficients  at  that  time.  ’ '  A  potential  drawback  is  that 
the  single-integral  inverse  does  not  represent  an  orthogonal  projection  under  the  energy 
metric.  However,  to  date  the  author  has  been  unable  distinguish  any  negative  impact. 


use  the  qualifier  quasi  because  this  does  not  remain  true  under  iteration.  That  is,  the  corre¬ 
sponding  exact  left  inverse  is  not  a  local  operator. 
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APPENDIX  A.  MORLET  WAVELETS 


The  generating  Morlet  wavelet  (of  norm  one)  tahes  the  form 


v>(t)  = 


TT 


1/4 


(A.!) 


There  are  two  parameters,  v  and  /3;  however,  in  general  only  the  ratio  fi/v  is  significant 
since  the  family 

Vit)  =  4=  ^  ^  (^.2) 

y'a  TT*/* 

obtained  by  scaling  only  depends  on  The  Fourier  transform  of  xl>{t)  is 

i?(u;)  =  tt'/"  ,  (^.3) 


and  the  energy  normalization  constants  and  Ci^  of  equations  (1.1)  and  (4.16)  are  given 

by 


and 


We  define  the  half-bandwidth,  Aw,  and  the  half-duration.  At,  to  be  the  point  at  which 
IV'I  and  [ipl  fall  to  1/e  of  their  peak  value;  i.e., 


Aw  =  \/2  ft 

(A. 6a) 

At  =  y/2/ft  . 

(AM) 

The  area  occupied  in 

pheise  space  is  thus 

(2  Aw)  (2  At)  =  8.0  271  . 

(A.7) 

Note  that  the  ratio 

=  If" 

At 

(A.8) 

» 


reflects  the  shape  of  the  region.  Smaller  /3  increases  frequency  resolution  while  decreasing 
time  resolution.  Since  (A. 8)  depends  on  /3^,  the  shape  is  fairly  sensitive  to  (j.  At  octave 
i  relations  (A.6)  become  Aw  =  /?  /  2’  and  At  =  2'  \/2  /  0.  Although  the  product  of 

these  two  remains  the  same,  their  ratio  is  Aw  /  At  =  Thus,  (A.8)  is  a  rather 

unsatisfying  description  of  /?.  A  more  invariant  and  intuitively  useful  description  of  0  is 
the  relative  frequency  resolution.  Since,  the  center  frequency  of  V’  is  v  /  2\  the  ratio  of 
the  bandwidth  to  the  central  frequency  is  {2\f20l2')  /  {v  j  2')  =  2\/2/3/n,  which  is 
independent  of  i. 


Discrete  Case 

For  convenience,  we  set  the  sample  rate  to  be  1.  Thus,  the  highpass  filter  g  is  de¬ 
termined  by  the  samples  ■>p{n).  Nyquist,  the  highest  frequency  representable,  is  w  =  tt. 
Hence,  in  order  to  avoid  aliasing  we  must  pick 

u  <  TT  —  \/20  .  (A. 9) 

For  0  =  0.5,  the  mzodmum  value  of  v  is  about  37r/4.  If  v  is  larger,  the  energy  will  alias, 
and  the  continuous  function  corresponding  to  g  will  not  be  the  Morlet  wavelet  That 
is,  the  DWT  will  be  the  wavelet  transform  with  respect  to  the  aliased  function.  Also,  in 
order  that  ip  be  admissible,  the  spectrum  must  be  essentially  0  for  w  <  0,  wliicl)  holds 
sufficiently  closely  if  [2]  -  [5] 

/3  <  ^  .  (/I.IO) 

To  achieve  this,  we  take  0  <  0.5  and  u  =  n  —  0.707. 

The  bandwidth  of  the  ujjper  half  of  the  spectrum  is  7r/2.  In  order  that  the  tran,sform 
retain  all  the  signal  energy,  and  hence  invert  in  a  numerically  stable  fashion,  the  highpass 
filters  must  cover  the  spectrum  from  7r/2  to  tt.  That  is,  we  must  have  enought  voices  to  do 
this.  If  7  =  2’/^^  where  L  is  the  number  of  voices,  then  voice  u  has  bandwidth  2  \/2  0  /  7". 
Thus,  the  total  bandwidth  covered  is 


Also, 


Hence, 


TDW  ^  2 

=  V?.0 


/?  f  1  +  ^  ’ 


^.2 

/ 


1-2-*/^^ 


1 


..L- 

I 


V  _  _L 

L— »oo  1  —  ln2 


TBW  « 


2L0  . 
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(A.ll) 

(A.12) 

(A.13) 
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To  cover  the  spectrum,  we  must  have  2L/3  >  7r/2;  i.e., 


L  > 


0.79 


1 

0 


(A14) 


which  is  generally  valid  for  13  <  0.25.  A  more  precise  formula  is  obtained  by  setting  (A.  11) 
greater  than  irJ2.  That  is,  /3  >  (tt  /  2  \/2)  (  1  —  2“^^^  )  which  yields, 


L  > 


In  2 

rn(l  -0.9i3)  ■ 


(A.15) 
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APPENDIX  B.  CONDITIONS  FOR  BOUNDED  TRANSFORMS 

We  derive  several  conditions  for  uniform  boundedness  (i.e.,  as  M  — f  oo)  of  the  DWT, 
the  DWT~^ ^  and  their  composition.  We  begin  with  necessary  conditions  expressed  as 
constraints  on  the  filters  at  a?  =  0,  and  finish  with  a  set  of  sufficient  conditions  (analogous  to 
(5.4))  for  boundedness  of  the  composition.  For  notational  clarity,  we  restrict  the  treatment 
to  a  single  voice. 

From  figure  2, 


w‘  =  (DV) »  n  *  1  = 


s’+’  =  [(D>/)  .  ]  s  , 


(5.1a) 


(5.16) 


with  the  obvious  interpretation  for  t  =  0  in  (B.la).  Using  (2.14),  we  obtain  the  z-transforms 
of  (B.l)  and  (5.6) 

j-i 

=  n  <7^(2'w)  s^(a))  (5.2a) 

j=0 


s"(w)  =  Y[  /*(2^^)Sz(w) 

7=0 

M-l  1-1  A/-1 

Sz(w)  =  X)  n  ':?z(2‘w)  w;(w)  +  ]][  7^(2>w)s^(w)  . 

i=0  7=0  7=0 

Composing  (B.2)  and  B.3),  we  have 

■M  — 1  J  — 1  t  — 1 

=  E  n  H  /,(2^u;)g,(2‘ca)j,,(2'ca) 

I-  1=0  7=0  7=0 

Af  — 1  M-]  •! 

I  TT  x  /07,  .^  TT  X  rni,  A  1., 

T  it  ii  J  z\^ 

7=0  7=0  J 

The  metric  induced  by  (5.2)  is 

M-i  , 

ii(w‘,s")ip  =  Yi  ■ 


{B.2b) 


(B.3) 


( n  A\ 


(5.5) 


It  follows  from  (B.2)  and  (B.3)  that  (5.5)  and  (5.7)  are  necessary  conditions  for  the  DWT 
and  DWT~^  to  be  bounded  with  respect  to  this  metric.  Also,  from  (B.4),  equation  (5.8) 


C*; 

W 

<*) 


I  • 
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is  a  necessary  condition  for  their  composition  to  be  bounded.  We  demonstrate  (5.5):  If 
the  DWT  is  bounded  (a,s  M  — »  oo)  there  must  exist  B  >  0)  such  that  for  all  M  and  s 

~  J  |s-(c.)|2  d<^  <  B  J  ls,(a;)|2  dLo  ^  B  l|sf  .  (B.6) 

Suppose  that  |/j(0)|  =  ■\/2(l  +  {3)  >  \/2.  Then,  there  exists  an  e  >  0  such  that  \f^(uj)\  > 
V2(l  +  a)  >  \/2  for  [wl  <  c.  Take  s  such  that  s*(w)  =  0  for  |tLi|  >  e/2'^.  Then 
|/^(2'w)l^  >  2(1  +  oY  for  i  <  M  and  Icj]  <  e/2’^  so  that  (B.2b)  implies  that 

^  /|<M“  >  i2"(l  +  af"  /  Is.HP 

=  (1+C.)2«||s.ll2  ,  (f?.7) 

which  contradicts  (B.6).  Thus,  for  boundedness,  we  must  have  1/^(0) |  <  \/2.  The  other 
proofs  utilize  <7^(0)  =  flFj,(0)  =  0. 

We  consider  the  composition  DWT~^  o  DWT  in  more  detail.  From  figure  2a,  we 
have  _ 

S'  =  (D'/)s^'  +  (D*S)w*  .  (B.8) 

Taking  the  z-transform  and  combining  it  with  that  of  fig;ire  2b  yields 

s;(u>)  =  J^{Tuj)sY\u2)  4  si,(2‘u.)w;(u.) 

=  7*(2'w)s-+1(u,)  +  ff,(2‘a7)ff,(2'a;)sl(c^)  .  (B.9) 

Lemma  Suppose  that 

+  C-^\9A<^)sY^)\  <  1  forallto  (5.10) 

where  C  >  1.  Then, 

|S‘»1  <  C|s'»|  .  (5.11) 

Proof  For  i  =  M,  we  have  s"(w)  =  s"(a;).  Assume  the  lemma  is  true  for  i  and  prove  it 

for  i  —  1. 


rsr‘(^)l  <  |/,(2’-'o;)S;(u;)l  +  C-Mff,(2’-Mf7.(2'-‘ca)CsrM 

<  |7,(2-ic.)Csi(a;)|  + 

=  l/,(2'-‘a:)/,(2-'a;)Csr(^)l  +  C-Mff,(2'-'a;)j,,(2’-'w)  Csri 

<  C'|sr‘(‘^-^)l  •  (5-12) 


A  C  >  1  satisfying  (B.IO)  exists  provided  \f  Y'^)f  z^^)\  —  ^ 


lgxfr)gz(‘^)l 


<  C  <  oo  for  all  u 


(5.13) 


(B.ll)  implies  that  ||sj|^  <  ||sll^.  If  g  is  complex  then,  a  forteriori,  ]|5es||^  <  lis|p. 

Thus,  (B.13)  is  a  sufficient  condition  for  DWT~^  o  DWT  to  be  a  bounded  transformation. 


APPENDIX  C.  NOTES  ON  DUALS 


The  ib*'*  column  of  a  matrix  represents  the  image  under  the  transformation  of  the 
basis  vector  in  the  domain.  Thus,  if  A  is  a  matrix  representing  the  forward  DWT,  its 
column  =  [Asj.^jt  i'’  the  image  of  [s]„  =  Stn-  But,  Wi^„  ~ 

Thus,  the  rows  of  A  correspond  to  the  Ijasis  vectors  ^^  „(*)•  Under  this  correspondence, 
we  have  - 

m 

In  other  words,  the  discrete  pseudo-inverse  ( A^  A)A^w  may  be  interpreted  as  an  expansion 
in  terms  of  dual  vectors  where  the  duals  V’i.nCO  given  by  (A^A)  V’i,n(') 

are  the  columns  of  Ah 
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(0.2) 


APPENDIX  D.  PLANCHEREL  RELATIONS  AND 
CONTINUOUS  INVERSION  FORMULAE 

Let  and  be  two  admissible  wavelet  functions  and  s’ (t)  and  s^{t)  two  signals 

in  L^.  Then  the  inner  product  of  the  corresponding  wavelet  transforms,  with  respect  to 
the  measure  da  db /a^,  is  equal  to  the  inner  product  of  the  two  signals  up  to  a  constant 
dependent  on  0’  and  That  is  [3], 

/oo  roo 

db  <s\  4>^l>  s'^  >  =  C^i^2<s\s‘^>  ,  (D.l) 

•OO  ^  J '-oo 

where  ^ 

J-oo 

If  we  abuse  the  requirement  that  be  (which  can  be  justified  if  we  consider  s“  to  be 
a  distribution)  and  take  s^(t)  =  S(t  —  r),  then  (D.l)  yields 

1  /In 

^\i)  =  r ^  4>K>4>Kit)  .  {D.3) 

J-oo  ^  J-oo 

If  xj)^  =  =  ?/>,  (D.3)  becomes  the  double  integration  inversion  formula 

1  da 

J^oo  ^  */— oo 

A 

where  and  =  <s,  xp^  >.  Relation  (D,4)  difTers  from  (1.5)  inasmuch  as  it 

includes  negative  values  of  a  (see  below).  Similarly,  if  ^^(t)  =  d{t)  we  obtain  the  single 

integration  formula.  It  is  clear  that  (D.2)  becomes  /  ^  |w|“’  tfw  =  Ci^  of  equation  (4.16), 
and  since 


(D.3)  becomes 


s  (t)  =  /  "7 

J-oo  i«r'‘ 


tp?  > 


=  r  Ti^  ■ 

Cirp  J-oo  m  ' 


(D.5) 


Formulae  for  a  >  0 

We  consider  formulae  for  positive  a.  Let  us  write  LHSa+  for  the  integral  of  the  left 
hand  side  of  (D.l),  but  with  a  restricted  to  [0,  oo].  The  derivation  in  [3]  yields 

/OO  ^  ^  ■—  /OO  J  _ _ ^ 

dw  s’(w)  s^(a;)  /  —  i/>’(au;)  (D.6) 

-oo  Jo  I®1 


•  % 
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p 


In  order  to  render  the  second  integral  independent  of  u>  through  the  substitution  a'  = 

we  need  to  have  x})^(a)rp'^{a)\a\~^ da  =  V’*(— a)V>2(— o)|a|“*da.  This  will  hold  true, 

for  example,  when  and  are  equal  and  real.  In  that  case,  the  result  is 

yoo  j  yOO  1 

/  — T  /  ,  ipl  >  <ipf ,  >  =  -C0<s’,s^>  with  tp  real  .  {D.7) 

Jo  “  J-CO  2 

We  observe  that  (D.7)  is  a  Plan chere- type  formula  stating  that  the  wavelet  transformation 
s  — >  (2/C^y^^W^  preserves  inner  products. 

Let 


<  s’,  > 


rH  tp\u)  p-  , 

(D.8a) 

V^’(w)  7/>2(a;)  ^  , 

-00  !^l 

(D.86) 

yOO  ^ 

=  /  dw  s’  s2  , 

Jo 

(D.9a) 

=  /  dw  s’  s2  . 

^  —00 

(DM) 

S^  >■*■  +  C^I0S<s’,S^> 

(D.IO) 

Then  (D.6)  may  be  rewritten 


If  s’  or  s^  is  analytic,  then  <  s’ ,  s'*  >  =0  and  <  s’ ,  s^  >  =  <  s’ ,  s’*  >■’"  so  that 
<5i^52>_  1 —  f  ^  f  db  <s\  rjJlxtp^l,  s'^  >  .  (-D-11) 

If  i/’’  or  ip'^  is  analytic  then  =  0  and  so 

r ^  r  rf/,  s^> .  (D.iz) 

Finally,  if  in  addition  to  the  analyticity  of  or  s’  and  s^  are  both  real,  we  have 
<  s’ ,  s^  >“  =  <  s’ ,  s2  >+  so  that  (D.12)  becomes 


<  s’,  s^  > 


=  iZe  I —  /  ^  /  d6  <s’,  V’t  >  <V'^^  >1  ■  (T>.i3) 

®  J  ~oo  J 


Taking  V>’  =  ip'^  =  ip  in  the  above  forrnulae,  as  well  as  s^(t)  =  6{t  —  r),  we  obtain  the 
inversion  formulae 
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db  V’jC*)  ^ 


(D.14) 


,  ,  2  /■”  <ia  /“ 

2  y*  oo  y«  oo 

s{t)  —  Re  I  —  I  db  '*>  analytic  tp  (JO.  15) 

Jq  a  J_^ 

1  /Ira 

s{t)  =  ~2  ^4*  analytic  s  .  (D.16) 

^ijf  Jo  ^  */'Oo 

To  get  the  single  integration  formulae,  we  take  —  S(t)  in  (D.6)  obtaining 

/oo  1  /OO  _  aOO  ^  — -  /oo  J  _ 

J  db  <s\  ip^>  s^(b)  =  J  (Lj  S^{uj)  s^{u<)  •  (^-17) 

Defining 

cu  -  io.iB) 

and  Cj“^  similarly,  we  have 

/OO  J  /OO 

Jo  J  v-?  >  =  C+^<s\s^>+  i  C-^<s\s‘^>-  .  (D,19) 

If  either  s'  or  s^  is  analytic,  this  yields 

1  da  _ 

<s',s^>  =  — ^  j  — -j^— ■  /  db  <  $^ ,  xp^  >  s^{b)  with  or  analytic  . 

^1  ti  Jo  1^1  J  —  00 


If  xp  is  analytic 


(D.20) 


1  da  _ 

<s',  s^>''’  =  —  /  -  -  ■ /  db  <  ,  xp^  >  s'^{b)  with  xp  analytic  ,  (D.21) 

Jo  |<ll  '  J-oo 

If  s*  and  are  real,  and  ^*(0;)  is  reed,  then  and  are  real  so  that 


LHS  +  LHS  =  (C+^  +  C,-^)  <s',  s2>  . 


(D.22) 


This  yields 


<  s',  Re 

Ci^ 


da  r 

Jo  \aY>^  Ic 


db  <  ,  xpl  >  3^ {b)  with  real  ,  xp  . 

(D.23) 


51 


The  corresponding  inversion  formulae  are 


s{t)  = 


<s,  rpi  >  real  s,  real  ^(w) 


(I>.24) 


s{t)  —  Re 


f  2  da  „  ) 

—  J  <s,  V*”  >  J  real  s,  analytic ‘4> 

1  da 

(t)  -  Tqr  /  tITT?  >  analytics  . 


(£>.25) 

(i:).26) 


Note  that  in  (D,25),  is  generally  complex.  Also,  there  seems  to  be  no  simple  single 
integration  formula  for  real  s{t)  with  real  ip(t). 


APPENDIX  E.  DERIVATION  OF  (5.19)  and(5.20) 


In  this  appendix,  it  is  convenient  as  well  as  enlightening  to  base  our  derivations  on  an 
alternative  implementation  to  that  found  in  figure  2;  one  which,  despite  appearances,  is 
more  directly  related  to  the  decimated  transform  of  figure  1.  Since  wj  =  Wi^o,  and  since  w* 
is  translation  invari2mt,  one  may  compute  the  undecimated  DWT,  wj,,  by  translating  the 
signal  by  n  and  then  computing  w,',o  (cf.  (2.9)).  This  may  be  efficiently  implemented  by 
the  structure  pictured  in  figure  5  [2].  This  diagram  is  most  easily  visualized  as  2^  copies 
of  figure  1 .  The  even  and  odd  branches  represent  a  splitting  of  the  time  series  into  its  even 
and  odd  samples.  Each  horizontal  path  through  the  diagram  is  identical  to  a  decimated 
DWT  filter  bank,  and  represents  a  time  advance  of  the  signal  by  ^  ■  2*ej  where  ti  is  1  for 
odd  branches  and  0  for  even  branches.  .4.  given  path  gets  repeated  every  2^  samples.  If 
we  index  the  paths  by  r,  we  see  that  as  r  goes  from  0  to  2'  —  1,  the  paths  correspond  to 
the  sets  of  wavelet  coefficients  Wi(r). 


Figure  5. 
transforms. 


Implementation  of  undecimated  wavelet  transform  by  muUiple.xing  decimated 


The  implementation  of  the  DWT~^  corresponding  to  figure  5,  is  obtained  by  simply 
inverting  its  data  flow.  In  this  direction,  the  outputs  of  the  filters  /  and  g  are  summed 
while  the  junctures  of  even  and  odd  branches  are  to  be  interpreted  as  multiplexing  the 
two  time  series.  This  multiplexing  can  be  represented  by  inserting  a  zero  in  between  the 
samples  of  each  channel,  delaying  the  odd  channel  by  one,  and  adding  the  two  channels. 
The  insertions  are  equivalent  to  applying  the  dilation  operator  D.  This  accounts  for  the 
delay  difference  of  one  sample  between  adjacent  channels  resulting,  ultimately,  in  a  delay  of 
r  samples  for  path  r  relative  to  path  0.  At  octave  i,  /*D  will  be  repeated  i  times  producing 
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a  contribution  at  time  m  of  the  form  [(/  +  U)*5  *  Wi(r)],„_r  from  path  r.  Adding  these 
contributions,  we  find 

2*^-1 

=  £  E  1(7 *D)'(J ♦■»,«)]„_,  +  E  l(7*D)"s“(r)],„-,  .  (£.1) 

i=0  r=0  r=0 

This  is  equation  (5.19).  To  prove  that  it  actually  computes  the  DWT~^ ^  we  must  demon¬ 
strate  its  equivalence  to  (5.6), 

Let  us  take  the  z-transform  of  (E.l).  We  first  observe  that 


=  53  =  S53^  ''Wi,„(r)2r  2’ «  , 

(£.2) 

k 

r=0  n  r=0  n 

which  yields 

2‘-l 

<(w)  =  53  • 

(£.3) 

Also,  one  has,  for  arbitrary  vectors  x  and  y, 

(y  *  Da:)z  =  yz(t^)  x^(2w)  . 

(£.4) 

Repeating  this  t  times. 

,  we  find 

i  — 1 

[(ytD)‘x].  =  ny,(2^a;)x,(2‘ai) 
j=o 

:-l 

(£.5) 

;=0 


Substituting  g  +  Wi(r)  for  x,  /  for  y,  and  summing  over  r,  we  obtain 

2'-l  2‘-l  i-I 

53  *  ^y(9  *  w.(r))],  =  ?z(2'w)  w,(r)^(2‘w) 

r=0  r=0 

1-1 

=  HiU^Uu,)  (D-s),  w;(w)  .  (E.Q) 

j  =  0 

Thi.s  is  the  z-transform  of  the  contribution  of  w'  to  s  (cf.,  (5.6)).  Summing  over  octaves  and 
including  the  DC  term,  we  obtain  equation  (5.6);  that  is,  the  algorithm  of  figure  2b.  This 
proves  (5.19).  To  complete  the  derivation  of  (5.20),  we  note  that  [(/Hxjni  =  'il>{rn~-n)xn- 

Replacing  V”  in  (5.18)  by  ^  ""  ’^)  9  (5.18)  by  *x  where 

[y']m  =  '</’'("r),  we  have 

[(/*D)*(i7^  ♦x)],„  Si  -  nyxn  =  Y^1p,^nXn  ■  {E.7) 

n  ^  ^  n 
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